Symbolic Combinatorics

3287 days ago by grady.thomas

# Grady Thomas, Ruoyu Zhang # MTHSC 865: Course Project # 06 December 2012 R.<z> = LazyPowerSeriesRing(QQ) #-------------------------------------------------------------------------- # The following defines a Combinatorial Class. Our implementation is only # concerned with obtaining a generating function from a given specification # so the class is essentially a container for the generating series. #-------------------------------------------------------------------------- class combinatorial_class: def __init__(self): self.generatingFunction = None def set_gf(self, coeff_list): #initialize generating functions of basic objects self.generatingFunction = R(coeff_list) def __add__(self, c): # Disjoint Union result_class = combinatorial_class() result_class.set_gf(self.generatingFunction + c.generatingFunction) return result_class sum = __add__ def __mul__(self, c): # Cartesian Product (Un-labelled Universe) & Labelled Product (Labelled Universe) result_class = combinatorial_class() result_class.set_gf(self.generatingFunction * c.generatingFunction) return result_class product = __mul__ #------------------------------------------------------------------------------------- # The following two classes are specific instances of a combinatorial class. These # classes are the basic building blocks of combinatorial objects. The neutral class # is a Combinatorial Class consisting of a single object of size 0. The atomic class # is a Combinatorial Class consisting of a single object of size 1. #------------------------------------------------------------------------------------- class neutral_class(combinatorial_class): def __init__(self): self.set_gf([1,0]) class atomic_class(combinatorial_class): def __init__(self): self.set_gf([0,1,0]) #----------------------------------------------------------------------------- # The following are a few methods we need for working with lazy power series. #------------------------------------------------------------------------------ def zk(n): # f(z) = z^k for i in range(n): yield 0 yield 1 yield 0 def pset_gen(f): # sum_{k=1}^{\infty} (-1)^(k-1)*(1/k)*f(z^k) n = 1 g = R([0]) while True: zn = R(zk(n)) s = (-1)^(n-1) g = s*(1/n)*f(zn) yield g n += 1 def mset_gen(f): # sum_{k=1}^{\infty} (1/k)*f(z^k) n = 1 g = R([0]) while True: zn = R(zk(n)) g = (1/n)*f(zn) yield g n += 1 #---------------------------------------------------------------------------- # Productions # A grammar specification of a combinatorial class is a statement of the form # C = rhs where C is an indicator of the class being named and rhs is one or # more of {E, Z, seq, pset, mset} where E and Z represent indicators for # the neutral and atomic classes, respectively. The following code handles # the translation of a specification to a generating function for the objects # constructable from that specification. #---------------------------------------------------------------------------- def seq(c, restriction='', n=0): # Sequences if restriction == '': sequences = R([1]) result_class = combinatorial_class() result_class.generatingFunction = sequences(c.generatingFunction) return result_class elif restriction == '=': # sequences having exactly n components result_class = combinatorial_class() f = c.generatingFunction for i in range(n-1): f = f*f result_class.generatingFunction = f return result_class elif restriction == '>=': # sequences having at least n components rclass = combinatorial_class() f = c.generatingFunction for i in range(n-1): f = f*f rclass.generatingFunction = f result_class = combinatorial_class() result_class = rclass*seq(c) return result_class def pset(c): # Power Set result_class = combinatorial_class() term = R.sum_generator(pset_gen(c.generatingFunction)) result_class.generatingFunction = term.exponential() return result_class def mset(c): # Multi-Set result_class = combinatorial_class() term = R.sum_generator(mset_gen(c.generatingFunction)) result_class.generatingFunction = term.exponential() return result_class 
       
# Some exampes of counting via specification. #------------------------------ # Words over the alphabet {a,b} #------------------------------ a = atomic_class() b = atomic_class() Words = seq(a+b) Words.generatingFunction.coefficients(20) 
       
\newcommand{\Bold}[1]{\mathbf{#1}}\left[1, 2, 4, 8, 16, 32, 64, 128, 256, 512, 1024, 2048, 4096, 8192, 16384, 32768, 65536, 131072, 262144, 524288\right]
\newcommand{\Bold}[1]{\mathbf{#1}}\left[1, 2, 4, 8, 16, 32, 64, 128, 256, 512, 1024, 2048, 4096, 8192, 16384, 32768, 65536, 131072, 262144, 524288\right]
#---------------------------------------------------- # Words over {a,b} matched by the regular expression # a((abb|babb(abb)+)+)* #---------------------------------------------------- p1 = a*seq(b, '=', 2) p2 = b*seq(p1, '>=', 2) rWords = a*seq(seq(p1 + p2, '>=', 1)) rWords.generatingFunction.coefficients(30) 
       
\newcommand{\Bold}[1]{\mathbf{#1}}\left[0, 1, 0, 0, 1, 0, 0, 2, 1, 0, 4, 5, 0, 8, 17, 2, 16, 49, 16, 32, 129, 78, 68, 321, 300, 172, 769, 1002, 536, 1801\right]
\newcommand{\Bold}[1]{\mathbf{#1}}\left[0, 1, 0, 0, 1, 0, 0, 2, 1, 0, 4, 5, 0, 8, 17, 2, 16, 49, 16, 32, 129, 78, 68, 321, 300, 172, 769, 1002, 536, 1801\right]
#---------------------------------------------- # The Natural numbers as a combinatorial class. #---------------------------------------------- Z = atomic_class() N = seq(Z, '>=', 1) N.generatingFunction.coefficients(15) 
       
\newcommand{\Bold}[1]{\mathbf{#1}}\left[0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1\right]
\newcommand{\Bold}[1]{\mathbf{#1}}\left[0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1\right]
#------------------------ # Compositions #------------------------ Compositions = pset(N) Compositions.generatingFunction.coefficients(15) 
       
\newcommand{\Bold}[1]{\mathbf{#1}}\left[1, 1, 1, 2, 2, 3, 4, 5, 6, 8, 10, 12, 15, 18, 22\right]
\newcommand{\Bold}[1]{\mathbf{#1}}\left[1, 1, 1, 2, 2, 3, 4, 5, 6, 8, 10, 12, 15, 18, 22\right]
#-------------------------------------------------------- # Compositions into at at least three components. #-------------------------------------------------------- Comps_3 = seq(N, '>=', 3) Comps_3.generatingFunction.coefficients(15) 
       
\newcommand{\Bold}[1]{\mathbf{#1}}\left[0, 0, 0, 0, 1, 5, 16, 42, 99, 219, 466, 968, 1981, 4017, 8100\right]
\newcommand{\Bold}[1]{\mathbf{#1}}\left[0, 0, 0, 0, 1, 5, 16, 42, 99, 219, 466, 968, 1981, 4017, 8100\right]
#-------------------------------------------------------- # Compositions into exactly three components. #-------------------------------------------------------- C3 = seq(N, '=', 3) C3.generatingFunction.coefficients(15) 
       
\newcommand{\Bold}[1]{\mathbf{#1}}\left[0, 0, 0, 0, 1, 4, 10, 20, 35, 56, 84, 120, 165, 220, 286\right]
\newcommand{\Bold}[1]{\mathbf{#1}}\left[0, 0, 0, 0, 1, 4, 10, 20, 35, 56, 84, 120, 165, 220, 286\right]
#------------------------------------------ # Partitions (Ordered Compositions) #------------------------------------------ Partitions = mset(N) Partitions.generatingFunction.coefficients(15) 
       
\newcommand{\Bold}[1]{\mathbf{#1}}\left[1, 1, 2, 3, 5, 7, 11, 15, 22, 30, 42, 56, 77, 101, 135\right]
\newcommand{\Bold}[1]{\mathbf{#1}}\left[1, 1, 2, 3, 5, 7, 11, 15, 22, 30, 42, 56, 77, 101, 135\right]
#----------------------------------------------------------------------------------------- # The implementation we have given currently only handles # iterative (non-recursive) specifications. A more thorough # implementation would be able to find generating functions for # objects which are recursively defined. # # In addition, classes may be 'marked' or 'weighted' in such a way that the coefficients in # the resulting power series are polynomials whose coefficients count # objects by both size and some other parameter. # # As a closing example we count (unlabelled) ordered trees by both total number of vertices and total # number of internal vertices using the Sage package for combinatorial species. #------------------------------------------------------------------------------------------- y = QQ['y'].gen() leaf = species.SingletonSpecies() internal_vertex = species.SingletonSpecies(weight=y) L = species.LinearOrderSpecies(min=1) T = species.CombinatorialSpecies() T.define(leaf+internal_vertex*L(T)) T.isotype_generating_series().coefficient(7) 
       
\newcommand{\Bold}[1]{\mathbf{#1}}y^{6} + 15 y^{5} + 50 y^{4} + 50 y^{3} + 15 y^{2} + y
\newcommand{\Bold}[1]{\mathbf{#1}}y^{6} + 15 y^{5} + 50 y^{4} + 50 y^{3} + 15 y^{2} + y
''' This tells us that, among all 132 ordered trees on 7 vertices: 1 tree has only a single internal vertex 15 trees have 2 internal vertices 50 trees have 3 internal vertices 50 trees have 4 internal vertices 15 trees have 5 internal vertices 1 tree has 6 internal vertices '''