Sage and Sage-Combinat demo, LACIM 2010, Montreal -- Sage Reference Manual v4.5 system:sage

Navigation

Sage and Sage-Combinat demo, LACIM 2010, Montreal

Basic calculations

{{{id=0| %hide pretty_print_default() /// }}}

Arithmetic:

{{{id=1| 1 + 1 /// }}} {{{id=2| factor(x^100 - 1) /// }}}

Symbolic calculations:

{{{id=3| var('y') f = sin(x) - cos(x*y) + 1 / (x^3+1) f /// }}} {{{id=4| f.integrate(x) /// }}}

Statistics:

{{{id=5| print r.summary(r.c(1,2,3,111,2,3,2,3,2,5,4)) /// }}}
sage: plot(sin(x), -pi, pi, fill = ‘axis’)

Taylor approximation:

{{{id=6| %hide var('x') @interact def g(f=sin(x), c=0, n=(1..30), xinterval=range_slider(-10, 10, 1, default=(-8,8), label="x-interval"), yinterval=range_slider(-50, 50, 1, default=(-3,3), label="y-interval")): x0 = c degree = n xmin,xmax = xinterval ymin,ymax = yinterval p = plot(f, xmin, xmax, thickness=4) dot = point((x0,f(x=x0)),pointsize=80,rgbcolor=(1,0,0)) ft = f.taylor(x,x0,degree) pt = plot(ft, xmin, xmax, color='red', thickness=2, fill=f) show(dot + p + pt, ymin=ymin, ymax=ymax, xmin=xmin, xmax=xmax) html('$f(x)\;=\;%s$'%latex(f)) html('$P_{%s}(x)\;=\;%s+R_{%s}(x)$'%(degree,latex(ft),degree)) /// }}}

Documentation

Type plot( and then press the TAB key:

{{{id=7| plot( /// }}}

Click on Help at the top of this page

There is also Live Documentation in which you can evaluate the examples.

Introspection:

{{{id=8| m = matrix([[1/2,1],[2,1]]); m /// }}} {{{id=9| m. /// }}}

Looking at the sources:

{{{id=10| m.determinant /// }}}

Algebraic constructions and categories

{{{id=11| Px.<x> = QQ[] Fx = Px.fraction_field() /// }}}
{{{id=12| for category in Fx.categories(): print category /// }}}
{{{id=13| g = sage.categories.category.category_graph() g.set_latex_options(format = "dot2tex") view(g, tightpage = True, viewer = "pdf") /// }}}

Cython: Python \longmapsto C

Here is a function that computes \sum_{k=0}^N k in pure Python:

{{{id=14| def mysum(N): s = int(0) for k in range(1,N): s += k return s time mysum(10^7) /// }}}

Let us compare this with the Cython version:

{{{id=15| %cython def mysum_cython(N): cdef long long s = 0 cdef int k for k in range(1,N): s += k return s time mysum_cython(10^7) /// }}}

A function to count the number of integer partitions with parts in a given set:

{{{id=16| def buying(coins, total): vlist = [ [0] * len(coins) for _ in range(total + 1) ] for i in range(total + 1): for j, coin in enumerate(coins): if j == 0: if i % coin == 0: vlist[i][j] = 1 else: k = 0 while i - k >= 0: vlist[i][j] += vlist[i - k][j - 1] k += coin return vlist[total][len(coins) - 1] [buying([1,2,5,10], i) for i in [1..20]] /// [1, 2, 2, 3, 4, 5, 6, 7, 8, 11, 12, 15, 16, 19, 22, 25, 28, 31, 34, 40] }}} {{{id=17| [1,3..10] /// [1, 3, 5, 7, 9] }}}

Let’s see how long it takes to find the number of partitions of 500 into odd parts:

{{{id=18| time buying([1,3..500], 500) /// 732986521245024 Time: CPU 3.05 s, Wall: 3.05 s }}}

Make two changes:

{{{id=19| %cython def cybuying(coins, total): cdef int i, j, k, coin vlist = [ [0] * len(coins) for _ in range(total + 1) ] for i in range(total + 1): for j, coin in enumerate(coins): if j == 0: if i % coin == 0: vlist[i][j] = 1 else: k = 0 while i - k >= 0: vlist[i][j] += vlist[i - k][j - 1] k += coin return vlist[total][len(coins) - 1] /// }}}

Surely two tiny changes in some Python code can’t make it much faster:

{{{id=20| time cybuying([1,3..500], 500) /// 732986521245024L Time: CPU 0.08 s, Wall: 0.09 s }}} {{{id=21| 3.05/0.08 /// 38.1250000000000 }}}

Sage-Combinat demonstration

{{{id=22| %hide pretty_print_default() /// }}}

Elementary combinatorics

Combinatorial objects

{{{id=23| p = Partition([3,3,2,1]) p /// }}} {{{id=24| p.pp() /// }}} {{{id=25| p.conjugate().pp() /// }}} {{{id=26| p.conjugate /// }}}
{{{id=27| s = Permutation([5,3,2,6,4,8,9,7,1]) s /// }}} {{{id=28| (p,q) = s.robinson_schensted() /// }}} {{{id=29| p.pp() /// 1 4 7 9 2 6 8 3 5 }}} {{{id=30| q.pp() /// 1 4 6 7 2 5 8 3 9 }}} {{{id=31| p.row_stabilizer() /// Permutation Group with generators [(), (7,9), (6,8), (4,7), (2,6), (1,4)] }}}

Enumerated sets (combinatorial classes)

{{{id=32| P5 = Partitions(5) P5 /// Partitions of the integer 5 }}} {{{id=33| P5.list() /// [[5], [4, 1], [3, 2], [3, 1, 1], [2, 2, 1], [2, 1, 1, 1], [1, 1, 1, 1, 1]] }}} {{{id=34| P5.cardinality() /// 7 }}} {{{id=35| Partitions(100000).cardinality() /// 27493510569775696512677516320986352688173429315980054758203125984302147328114964173055050741660736621590157844774296248940493063070200461792764493033510116079342457190155718943509725312466108452006369558934464248716828789832182345009262853831404597021307130674510624419227311238999702284408609370935531629697851569569892196108480158600569421098519 }}} {{{id=36| Permutations(20).random_element() /// [15, 6, 8, 14, 17, 16, 4, 7, 11, 3, 10, 5, 19, 9, 12, 2, 20, 18, 1, 13] }}} {{{id=37| Compositions(10).unrank(100) # TODO: non stupid algorithm /// [1, 1, 3, 1, 2, 1, 1] }}} {{{id=38| for p in StandardTableaux([3,2]): print "-----------------------------" p.pp() /// ----------------------------- 1 3 5 2 4 ----------------------------- 1 2 5 3 4 ----------------------------- 1 3 4 2 5 ----------------------------- 1 2 4 3 5 ----------------------------- 1 2 3 4 5 }}}

Trees

ToDo

Summary:

  • Every mathematical object (element, set, category, ...) is modeled by a Python object</li>
  • All combinatorial classes share a uniform interface</li>

Constructions

{{{id=39| C = DisjointUnionEnumeratedSets( [ Compositions(4), Permutations(3)] ) C /// Union of Family (Compositions of 4, Standard permutations of 3) }}} {{{id=40| C.cardinality() /// 14 }}} {{{id=41| C.list() /// [[1, 1, 1, 1], [1, 1, 2], [1, 2, 1], [1, 3], [2, 1, 1], [2, 2], [3, 1], [4], [1, 2, 3], [1, 3, 2], [2, 1, 3], [2, 3, 1], [3, 1, 2], [3, 2, 1]] }}}
{{{id=42| C = CartesianProduct(Compositions(8), Permutations(20)) C /// Cartesian product of Compositions of 8, Standard permutations of 20 }}} {{{id=43| C.cardinality() /// 311411457046609920000 }}} {{{id=44| C.random_element() # todo: broken /// }}}
{{{id=45| F = Family(NonNegativeIntegers(), Permutations) F /// Lazy family (Permutations(i))_{i in Set of non negative integers} }}} {{{id=46| F[1000] /// Standard permutations of 1000 }}} {{{id=47| U = DisjointUnionEnumeratedSets(F) U.cardinality() /// +Infinity }}} {{{id=48| p = iter(U) for i in range(12): print p.next() /// [] [1] [1, 2] [2, 1] [1, 2, 3] [1, 3, 2] [2, 1, 3] [2, 3, 1] [3, 1, 2] [3, 2, 1] [1, 2, 3, 4] [1, 2, 4, 3] }}} {{{id=49| for p in U: print p /// [] [1] [1, 2] [2, 1] [1, 2, 3] [1, 3, 2] [2, 1, 3] [2, 3, 1] [3, 1, 2] }}}

Summary:

  • Basic combinatorial classes + constructions give a flexible toolbox
  • This is made possible by uniform interfaces
  • Lazy algorithms and data structures for large / infinite sets (iterators, ...)

Enumeration kernels

Integer lists:

{{{id=50| IntegerVectors(10, 3, min_part = 2, max_part = 5, inner = [2, 4, 2]).list() /// [[4, 4, 2], [3, 5, 2], [3, 4, 3], [2, 5, 3], [2, 4, 4]] }}} {{{id=51| Compositions(5, max_part = 3, min_length = 2, max_length = 3).list() /// [[1, 1, 3], [1, 2, 2], [1, 3, 1], [2, 1, 2], [2, 2, 1], [2, 3], [3, 1, 1], [3, 2]] }}} {{{id=52| Partitions(5, max_slope = -1).list() /// [[5], [4, 1], [3, 2]] }}} {{{id=53| IntegerListsLex(10, length=3, min_part = 2, max_part = 5, floor = [2, 4, 2]).list() /// [[4, 4, 2], [3, 5, 2], [3, 4, 3], [2, 5, 3], [2, 4, 4]] }}} {{{id=54| IntegerListsLex(5, min_part = 1, max_part = 3, min_length = 2, max_length = 3).list() /// [[3, 2], [3, 1, 1], [2, 3], [2, 2, 1], [2, 1, 2], [1, 3, 1], [1, 2, 2], [1, 1, 3]] }}} {{{id=55| IntegerListsLex(5, min_part = 1, max_slope = -1).list() /// [[5], [4, 1], [3, 2]] }}} {{{id=56| c = Compositions(5)[1] c /// [1, 1, 1, 2] }}} {{{id=57| c = IntegerListsLex(5, min_part = 1)[1] /// }}}

Species / decomposable classes

{{{id=58| from sage.combinat.species.library import * o = var("o") /// }}}

Fibonacci words:

{{{id=59| Eps = EmptySetSpecies() Z0 = SingletonSpecies() Z1 = Eps*SingletonSpecies() FW = CombinatorialSpecies() FW.define(Eps + Z0*FW + Z1*Eps + Z1*Z0*FW) FW /// }}} {{{id=60| L = FW.isotype_generating_series().coefficients(15) L /// }}} {{{id=61| sloane_find(L) /// Searching Sloane's online database... [[45, 'Fibonacci numbers: F(n) = F(n-1) + F(n-2), F(0) = 0, F(1) = 1, F(2) = 1, ...', [0, 1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, 610, 987, 1597, 2584, 4181, 6765, 10946, 17711, 28657, 46368, 75025, 121393, 196418, 317811, 514229, 832040, 1346269, 2178309, 3524578, 5702887, 9227465, 14930352, 24157817, 39088169]], [24595, 'a(n) = s(1)t(n) + s(2)t(n-1) + ... + s(k)t(n+1-k), where k = [ (n+1)/2 ], s = (F(2), F(3), ...), t = A023533.', [1, 0, 0, 1, 2, 3, 5, 0, 0, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, 610, 987, 1598, 2586, 4184, 6770, 10954, 13, 21, 34, 55, 89, 144, 233, 377, 610, 987, 1597, 2584, 4181, 6765, 10946, 17711, 28658, 46370, 75028, 121398, 196426]], [25109, 'a(n) = s(1)t(n) + s(2)t(n-1) + ... + s(k)t(n-k+1), where k = [ n/2 ], s = (F(2), F(3), F(4), ...), t = A023533.', [0, 0, 1, 2, 3, 0, 0, 0, 1, 2, 3, 5, 8, 13, 21, 34, 55, 0, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, 610, 987, 1598, 2586, 4181, 6770, 8, 13, 21, 34, 55, 89, 144, 233, 377, 610, 987, 1597, 2584, 4181, 6765, 10946, 17711, 28658, 46370, 75028, 121398, 196426, 317824, 514250]], [132636, 'Fib(n) mod n^3.', [0, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, 610, 987, 1597, 2584, 4181, 6765, 1685, 7063, 4323, 4896, 12525, 15937, 19271, 10483, 2060, 22040, 5674, 15621, 2752, 3807, 9340, 432, 46989, 19305, 11932, 62155, 31899, 12088, 22273, 3677, 32420]], [132916, 'a(0)=0; a(1)=1; a(n) = Sum a(n-k), k= 1 ... [n^(1/3)] for n>=2.', [0, 1, 1, 1, 1, 1, 1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, 610, 987, 1597, 2584, 4181, 6765, 10946, 21892, 39603, 72441, 133936, 245980, 452357, 832273, 1530610, 2815240, 5178123, 9523973, 17517336, 32219432, 59260741, 108997509, 200477682]], [147316, 'A000045 Fibonacci mirror sequence Binet: f(n)=(1/5)*2^(-n) ((5 - 2 *Sqrt[5]) (1 + Sqrt[5])^n + (1 - Sqrt[5])^n(5 + 2 * Sqrt[5])).', [1597, -987, 610, -377, 233, -144, 89, -55, 34, -21, 13, -8, 5, -3, 2, -1, 1, 0, 1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, 610, 987, 1597]], [39834, 'a(n+2)=-a(n+1)+a(n) (signed Fibonacci numbers); or Fibonacci numbers (A000045) extended to negative indices.', [1, 1, 0, 1, -1, 2, -3, 5, -8, 13, -21, 34, -55, 89, -144, 233, -377, 610, -987, 1597, -2584, 4181, -6765, 10946, -17711, 28657, -46368, 75025, -121393, 196418, -317811, 514229, -832040, 1346269, -2178309, 3524578, -5702887, 9227465, -14930352, 24157817]], [152163, 'a(n)=a(n-1)+a(n-2), n>1 ; a(0)=1, a(1)=-1 .', [1, -1, 0, -1, -1, -2, -3, -5, -8, -13, -21, -34, -55, -89, -144, -233, -377, -610, -987, -1597, -2584, -4181, -6765, -10946, -17711, -28657, -46368, -75025, -121393, -196418, -317811, -514229, -832040, -1346269, -2178309, -3524578, -5702887]]] }}} {{{id=62| FW3 = FW.isotypes([o]*4); FW3 /// Isomorphism types for Combinatorial species with labels [o, o, o, o] }}} {{{id=63| FW3.list() /// [o*(o*(o*(o*{}))), o*(o*(o*(({}*o)*{}))), o*(o*((({}*o)*o)*{})), o*((({}*o)*o)*(o*{})), o*((({}*o)*o)*(({}*o)*{})), (({}*o)*o)*(o*(o*{})), (({}*o)*o)*(o*(({}*o)*{})), (({}*o)*o)*((({}*o)*o)*{})] }}}

Binary trees:

{{{id=64| BT = CombinatorialSpecies() Leaf = SingletonSpecies() BT.define(Leaf+(BT*BT)) BT5 = BT.isotypes([o]*5) /// }}} {{{id=65| BT5.list() /// [o*(o*(o*(o*o))), o*(o*((o*o)*o)), o*((o*o)*(o*o)), o*((o*(o*o))*o), o*(((o*o)*o)*o), (o*o)*(o*(o*o)), (o*o)*((o*o)*o), (o*(o*o))*(o*o), ((o*o)*o)*(o*o), (o*(o*(o*o)))*o, (o*((o*o)*o))*o, ((o*o)*(o*o))*o, ((o*(o*o))*o)*o, (((o*o)*o)*o)*o] }}} {{{id=66| %hide def pbt_to_coordinates(t): e = {} queue = [t] while queue: z = queue.pop() if not isinstance(z[0], int): e[z[1]._labels[0]-1] = z queue.extend(z) coord = [(len(e[i][0]._labels) * len(e[i][1]._labels)) for i in range(len(e))] return sage.geometry.polyhedra.Polytopes.project_1(coord) K4 = Polyhedron(vertices=[pbt_to_coordinates(t) for t in BT.isotypes(range(5))]) K4.show(fill=True).show(frame=False) /// }}}

Juggling automaton:

{{{id=67| F = SingletonSpecies() state_labels = map(tuple, Permutations([0,0,1,1,1]).list()) states = dict((i, CombinatorialSpecies()) for i in state_labels) def successors(state): newstate = state[1:]+(0,) if state[0] == 0: return [(0, newstate)] return [(i+1, newstate[0:i] + (1,) + newstate[i+1:]) for i in range(0, len(state)) if newstate[i] == 0] for state in state_labels: states[state].define( sum( [states[target]*F for (height, target) in successors(state)], Eps )) states /// {(1, 1, 0, 1, 0): Combinatorial species, (0, 1, 1, 0, 1): Combinatorial species, (1, 1, 1, 0, 0): Combinatorial species, (1, 0, 1, 0, 1): Combinatorial species, (0, 1, 0, 1, 1): Combinatorial species, (1, 0, 0, 1, 1): Combinatorial species, (0, 1, 1, 1, 0): Combinatorial species, (1, 0, 1, 1, 0): Combinatorial species, (0, 0, 1, 1, 1): Combinatorial species, (1, 1, 0, 0, 1): Combinatorial species} }}} {{{id=68| states[(1,1,1,0,0)].isotypes([o]*5).cardinality() /// 165 }}}

Lattice points of polytopes

{{{id=69| A=random_matrix(ZZ,3,6,x=7) L=LatticePolytope(A) L.plot3d() /// }}} {{{id=70| L.npoints() # should be cardinality! /// 28 }}}

This example used PALP and J-mol

Graphs up to an isomorphism

{{{id=71| show(graphs(5, lambda G: G.size() <= 4)) /// }}}

Words

An infinite periodic word:

{{{id=72| p = Word([0,1,1,0,1,0,1]) ^ Infinity p /// word: 0110101011010101101010110101011010101101... }}}

The Fibonacci word:

{{{id=73| f = words.FibonacciWord() f /// word: 0100101001001010010100100101001001010010... }}}

The Thue-Morse word:

{{{id=74| t = words.ThueMorseWord() t /// word: 0110100110010110100101100110100110010110... }}}

A random word over the alphabet [0, 1] of length 1000:

{{{id=75| r = words.RandomWord(1000,2) r /// word: 0010101011101100110000000110000111011100... }}}

The fixed point of a morphism:

{{{id=76| m = WordMorphism('a->acabb,b->bcacacbb,c->baba') w = m.fixed_point('a') /// }}} {{{id=77| w /// word: acabbbabaacabbbcacacbbbcacacbbbcacacbbac... }}}

Their prefixes of length 1000:

{{{id=78| pp = p[:1000] ff = f[:1000] tt = t[:1000] ww = w[:1000] /// }}}

A comparison of their complexity function:

{{{id=79| A = list_plot([pp.number_of_factors(i) for i in range(100)], color='red') B = list_plot([ff.number_of_factors(i) for i in range(100)], color='blue') C = list_plot([tt.number_of_factors(i) for i in range(100)], color='green') D = list_plot([r.number_of_factors(i) for i in range(100)], color='black') E = list_plot([ww.number_of_factors(i) for i in range(100)], color='orange') A + B + C + D + E /// }}}

Construction of a permutation and builds its associated Rauzy diagram:

{{{id=80| p = iet.Permutation('a b c d', 'd c b a') p /// a b c d d c b a }}} {{{id=81| r = p.rauzy_diagram() r /// Rauzy diagram with 7 permutations }}} {{{id=82| r.graph().plot() /// }}}

Let us now construct a self-similar interval exchange transformation associated to a loop in the Rauzy diagram:

{{{id=83| g0 = r.path(p, 't', 't', 'b', 't') g1 = r.path(p, 'b', 'b', 't', 'b') g = g0 + g1 m = g.matrix() v = m.eigenvectors_right()[-1][1][0] T = iet.IntervalExchangeTransformation(p, v) /// }}}

We can plot it and all its power:

{{{id=84| T.plot() (T*T).plot() (T*T*T).plot() /// }}}

Check the self similarity of T:

{{{id=85| T.rauzy_diagram(iterations=8).normalize(T.length()) == T /// True }}}

And get the symbolic coding of 0 using the substitution associated to the path:

{{{id=86| s = g.orbit_substitution() s.fixed_point('a') /// word: adbdadcdadbdbdadcdadbdadcdadccdadcdadbda... }}}

Predefined algebraic structures

Root systems, Coxeter groups, ...

{{{id=87| L = RootSystem(['A',2,1]).weight_space() L.plot(size=[[-1..1],[-1..1]],alcovewalks=[[0,2,0,1,2,1,2,0,2,1]]) /// }}} {{{id=88| CartanType.samples() /// [['A', 1], ['A', 5], ['B', 1], ['B', 5], ['C', 1], ['C', 5], ['D', 2], ['D', 3], ['D', 5], ['E', 6], ['E', 7], ['E', 8], ['F', 4], ['G', 2], ['I', 5], ['H', 3], ['H', 4], ['A', 1, 1], ['A', 5, 1], ['B', 1, 1], ['B', 5, 1], ['C', 1, 1], ['C', 5, 1], ['D', 3, 1], ['D', 5, 1], ['E', 6, 1], ['E', 7, 1], ['E', 8, 1], ['F', 4, 1], ['G', 2, 1], ['B', 5, 1]^*, ['C', 4, 1]^*, ['F', 4, 1]^*, ['G', 2, 1]^*, ['BC', 1, 2], ['BC', 5, 2]] }}} {{{id=89| T = CartanType(["E", 8, 1]) print T.dynkin_diagram() /// O 2 | | O---O---O---O---O---O---O---O 1 3 4 5 6 7 8 0 E8~ }}} {{{id=90| T.is_simply_laced(), T.is_finite(), T.is_crystalographic() /// (True, False, True) }}} {{{id=91| W = WeylGroup(["B", 3]) W /// Weyl Group of type ['B', 3] (as a matrix group acting on the ambient space) }}} {{{id=92| W.cayley_graph(side = "left").plot3d(color_by_label = True) /// }}} {{{id=93| print W.character_table() # Thanks GAP! /// CT1 2 4 4 3 3 4 3 1 1 3 4 3 1 . . . . . 1 1 . 1 1a 2a 2b 4a 2c 2d 6a 3a 4b 2e X.1 1 1 1 1 1 1 1 1 1 1 X.2 1 1 1 -1 -1 -1 -1 1 1 -1 X.3 1 1 -1 -1 1 -1 1 1 -1 1 X.4 1 1 -1 1 -1 1 -1 1 -1 -1 X.5 2 2 . . -2 . 1 -1 . -2 X.6 2 2 . . 2 . -1 -1 . 2 X.7 3 -1 1 1 1 -1 . . -1 -3 X.8 3 -1 -1 -1 1 1 . . 1 -3 X.9 3 -1 -1 1 -1 -1 . . 1 3 X.10 3 -1 1 -1 -1 1 . . -1 3 }}} {{{id=94| rho = SymmetricGroupRepresentation([3, 2], "orthogonal"); rho /// Orthogonal representation of the symmetric group corresponding to [3, 2] }}} {{{id=95| rho([1, 3, 2, 4, 5]) /// 1 & 0 & 0 & 0 & 0 \\ 0 & -\frac{1}{2} & \frac{1}{2} \, \sqrt{3} & 0 & 0 \\ 0 & \frac{1}{2} \, \sqrt{3} & \frac{1}{2} & 0 & 0 \\ 0 & 0 & 0 & -\frac{1}{2} & \frac{1}{2} \, \sqrt{3} \\ 0 & 0 & 0 & \frac{1}{2} \, \sqrt{3} & \frac{1}{2} }}}

Affine Weyl groups:

{{{id=96| W = WeylGroup(["C", 3, 1]) W /// Weyl Group of type ['C', 3, 1] (as a matrix group acting on the root space) }}} {{{id=97| W.category() /// Category of affine weyl groups }}} {{{id=98| W.an_element() /// [-1 1 0 0] [ 0 1 0 0] [ 0 0 1 0] [ 0 0 0 1] }}} {{{id=99| W.from_reduced_word([1,2,3,0,3,0,3,2,1,3,3,2]).stanley_symmetric_function() /// 256*m[1, 1, 1, 1, 1, 1] + 128*m[2, 1, 1, 1, 1] + 64*m[2, 2, 1, 1] + 32*m[2, 2, 2] + 48*m[3, 1, 1, 1] + 24*m[3, 2, 1] + 8*m[3, 3] + 16*m[4, 1, 1] + 8*m[4, 2] + 4*m[5, 1] }}}

Symmetric functions

Classical basis:

{{{id=100| Sym = SymmetricFunctions(QQ) Sym /// Symmetric Functions over Rational Field }}} {{{id=101| s = Sym.schur() h = Sym.complete() e = Sym.elementary() m = Sym.monomial() p = Sym.powersum() /// }}} {{{id=102| m(( ( h[2,1] * ( 1 + 3 * p[2,1]) ) + s[2](s[3]))) /// }}}

Jack polynomials:

{{{id=103| Sym = SymmetricFunctions(QQ['t']) Jack = Sym.jack_polynomials() # todo: not implemented P = Jack.P(); J = Jack.J(); Q = Jack.Q() # todo: not implemented J(P[2,1]) # todo: not implemented /// AttributeError: 'SymmetricFunctions_with_category' object has no attribute 'jack_polynomials' }}}

Macdonald polynomials:

{{{id=104| J = MacdonaldPolynomialsJ(QQ) P = MacdonaldPolynomialsP(QQ) Q = MacdonaldPolynomialsQ(QQ) J /// Macdonald polynomials in the J basis over Fraction Field of Multivariate Polynomial Ring in q, t over Rational Field }}} {{{id=105| f = P(J[2,2] + 3 * Q[3,1]) f /// (q^2*t^6-q^2*t^5-q^2*t^4-q*t^5+q^2*t^3+2*q*t^3+t^3-q*t-t^2-t+1)*McdP[2, 2] + ((3*q^3*t^5-6*q^3*t^4+3*q^3*t^3-3*q^2*t^4+6*q^2*t^3-3*q^2*t^2-3*q*t^3+6*q*t^2-3*q*t+3*t^2-6*t+3)/(q^7*t-2*q^6*t+2*q^4*t-q^4-q^3*t+2*q^3-2*q+1))*McdP[3, 1] }}} {{{id=106| f /// }}} {{{id=107| Sym = SymmetricFunctions(J.base_ring()) s = Sym.s() s(f) /// }}}

Semigroups

{{{id=108| Cat = FiniteSemigroups() Cat /// Category of finite semigroups }}} {{{id=109| Cat.category_graph().plot(layout = "acyclic") /// }}} {{{id=110| S = Cat.example(alphabet = ('a','b','c')) S /// An example of finite semi-group: the left regular band generated by ('a', 'b', 'c') }}} {{{id=111| S.cayley_graph(side = "left", simple = True).plot() /// }}} {{{id=112| S.j_transversal_of_idempotents() /// ['cab', 'ca', 'ab', 'cb', 'a', 'c', 'b'] }}} {{{id=113| S /// }}}

Hopf algebras

{{{id=114| Cat = HopfAlgebrasWithBasis(QQ); Cat /// Category of hopf algebras with basis over Rational Field }}} {{{id=115| g = Cat.category_graph() g.set_latex_options(format = "dot2tex") view(g, tightpage = True, viewer = "pdf") /// }}} {{{id=116| Cat H = Cat.example() H /// The Hopf algebra of the Dihedral group of order 6 as a permutation group over Rational Field }}} {{{id=117| H /// }}}

A real life example

{{{id=118| def path_to_line(path, grid=True): vert = lambda x, y: circle((x, y), .05, rgbcolor=(0, 0, 1), fill=True) pline = [(0,0)] vertices = [vert(0,0)] h = 0 maxh = h ln = len(path) for x, y in enumerate(path): h += y if h > maxh: maxh = h pline += [(x+1, h)] vertices += [vert(x+1, h)] plotted_path = line(pline) + sum(vertices) if grid: gridline = lambda a, b, c, d: line([(a, b), (c,d)], rgbcolor=(0,) * 3, linestyle='--', alpha=.25) # vertical gridlines grid = [gridline(x, 0, x, maxh) for x in [1..ln]] # horizontal gridlines for y in [1..maxh]: grid += [gridline(0, y, ln, y)] plotted_path += sum(grid) plotted_path.set_aspect_ratio(1) return plotted_path from sage.combinat.backtrack import GenericBacktracker class LukPaths(GenericBacktracker): def __init__(self, n, k=1): GenericBacktracker.__init__(self, [], (0, 0)) self._n = n self._k = k if n < 0 or k < 1 or n % (k+1) != 0: def jane_stop_this_crazy_thing(*args): raise StopIteration self._rec = jane_stop_this_crazy_thing def _rec(self, path, state): ln, ht = state if ln < self._n: # if we're high enough, we can yield a path with a # k-downstep at the end if ht >= self._k: yield path + [-self._k], (ln + 1, ht - self._k), False # if the path isn't too high, it can also take an upstep if ht / self._k < (self._n - ln): yield path + [1], (ln + 1, ht + 1), False else: # if length is n, set state to None so we stop trying to # make new paths, and yield what we've got yield path, None, True plots = [path_to_line(p) for p in LukPaths(12, 2)] ga = graphics_array(plots, ceil(len(plots)/6), 6) ga.show(figsize=[9,7]) /// }}}

Real life example: Parallel testing of a conjecture on J-Trivial monoids using MuPAD

{{{id=119| %hide pretty_print_default() /// }}}
{{{id=120| from sage.combinat.j_trivial_monoids import * def pij(j): return lambda i: i if i != j+1 else j pi2 = pij(2) /// }}}
{{{id=121| pi2(1), pi2(2), pi2(3), pi2(4) /// (1, 2, 2, 4) }}}
{{{id=122| class NDPFMonoid(AutomaticMonoid): def __init__(self, n): ambient_monoid = DiscreteFunctions(range(1,n+1), action="right") pi = Family(range(1, n), lambda j: ambient_monoid(pij(j))) AutomaticMonoid.__init__(self, pi, one = ambient_monoid.one(), category = (SubFiniteMonoidsOfFunctions(), FiniteJTrivialMonoids())) Mon = NDPFMonoid(3) Mon.cardinality() /// 5 }}}
{{{id=123| Mon.list() /// [[], [1], [2], [1, 2], [2, 1]] }}}
{{{id=124| [ NDPFMonoid(n).cardinality() for n in range(6)] /// [1, 1, 2, 5, 14, 42] }}}
{{{id=125| MuMon = mupad(Mon); MuMon /// / +- -+ \ | | 0, 1, 2, 3, 4 | | | | | | | | 1, 1, 4, 4, 4 | | | | | | Dom::MMonoid| | 2, 3, 2, 3, 4 | | | | | | | | 3, 3, 4, 4, 4 | | | | | | | | 4, 4, 4, 4, 4 | | \ +- -+ / }}} {{{id=126| MuMon.count() /// 5 }}} {{{id=127| MuAlg = mupad.Dom.MonoidAlgebra(MuMon); MuAlg /// }}} {{{id=128| MuCMat = MuAlg.cartanInvariantsMatrix(); MuCMat /// +- -+ | 1, 0, 0, 0 | | | | 0, 1, 1, 0 | | | | 0, 0, 1, 0 | | | | 0, 0, 0, 1 | +- -+ }}} {{{id=129| MuCMat.sage() /// [1 0 0 0] [0 1 1 0] [0 0 1 0] [0 0 0 1] }}}
{{{id=130| M4 = NDPFMonoid(4) var('q') /// q }}} {{{id=131| cartconj = M4.cartan_matrix(q); cartconj /// [ 1 0 0 0 0 0 0 0] [ 0 1 q q^2 0 0 0 0] [ 0 0 1 q 0 0 0 0] [ 0 0 0 1 0 0 0 0] [ 0 0 0 0 1 0 q 0] [ 0 0 0 0 q 1 q^2 0] [ 0 0 0 0 0 0 1 0] [ 0 0 0 0 0 0 0 1] }}} {{{id=132| cart = M4.cartan_matrix_mupad(q); cart /// [ 1 0 0 0 0 0 0 0] [ 0 1 0 0 q 0 0 0] [ 0 q 1 0 q^2 0 0 0] [ 0 0 0 1 0 q^2 q 0] [ 0 0 0 0 1 0 0 0] [ 0 0 0 0 0 1 0 0] [ 0 0 0 0 0 q 1 0] [ 0 0 0 0 0 0 0 1] }}} {{{id=133| def is_isomorphic_matrices(m1, m2): coeffs1 = set([ c for row in m1 for c in row ]) coeffs2 = set([ c for row in m2 for c in row ]) if coeffs1 != coeffs2: return False f = sage.combinat.ranker.rank_from_list(sorted(coeffs1)) def graph(m): m = matrix([[f(m[i,j]) for j in range(m.ncols()) ] for i in range(m.nrows())]) return DiGraph(m, multiple_edges = True) return graph(m1).is_isomorphic(graph(m2)) /// }}} {{{id=134| is_isomorphic_matrices(cart, cartconj) /// True }}} {{{id=135| P4 = Posets(4); P4 /// Posets containing 4 vertices }}} {{{id=136| P4.cardinality() /// 16 }}} {{{id=137| Pos = P4[9]; Pos.cover_relations() /// [[0, 2], [1, 2], [2, 3]] }}} {{{id=138| #Pos.plot() /// }}} {{{id=139| MP = NDPFMonoidPoset(Pos); MP /// NDPF monoid of Poset ([[0, 2], [1, 2], [2, 3]]) }}} {{{id=140| is_isomorphic_matrices(MP.cartan_matrix(q), MP.cartan_matrix_mupad(q)) /// True }}} {{{id=141| @parallel() def check_conj_parallel(Pos): MP = NDPFMonoidPoset(Pos) return is_isomorphic_matrices(MP.cartan_matrix(q), MP.cartan_matrix_mupad(q)) /// }}} {{{id=142| for (((poset,), _), res) in check_conj_parallel(Posets(3).list()): print poset.cover_relations(), res /// }}} {{{id=143| all(x[1] for x in check_conj_parallel(Posets(4).list())) /// True }}}

Table Of Contents

This Page

Navigation

© Copyright 2005--2010, The Sage Development Team. Created using Sphinx 0.6.3.