Tutorial: Implementing Algebraic Structures

This tutorial will discuss five concepts:

At the end of this tutorial, the reader should be able to reimplement by himself the example of algebra with several realizations:

{{{id=0| Sets().WithRealizations() /// }}}

Namely, we consider an algebra $A(S)$ whose basis is indexed by the subsets $s$ of a given set $S$. $A(S)$ is endowed with three natural basis: F, In, Out; in the first basis, the product is given by the union of the indexing sets. The In basis and Out basis are defined respectively by:

Unknown directive type "math".

.. math::

    In_s  = \sum_{ t\subset s}  F_t \qquad
    F_s   = \sum_{ t\supset s}  Out_t


Each such basis gives a realization of $A$, where the elements are represented by their expansion in this basis. In the running exercises we will progressively implement this algebra and its three realizations, with coercions and mixed arithmetic between them.

This tutorial heavily depends on the `tutorial-using-free-modules`_.

Subclassing free modules and including category information

As a warm-up, we implement the group algebra of the additive group $ZZ/5ZZ$. Of course this is solely for pedagogical purposes; group algebras are already implemented (see ZMod(5).algebra(ZZ)). Recall that a fully functional $ZZ$-module over this group can be created with the simple command:

{{{id=1| A = CombinatorialFreeModule(ZZ, Zmod(5), prefix='a') /// }}}

We reproduce the same, but by deriving a subclass of :class:`CombinatorialFreeModule`:

{{{id=2| class MyCyclicGroupModule(CombinatorialFreeModule): """An absolutely minimal implementation of a module whose basis is a cyclic group""" def __init__(self, R, n, *args, **kwargs): CombinatorialFreeModule.__init__(self, R, Zmod(n), *args, **kwargs) /// }}} {{{id=3| A = MyCyclicGroupModule(QQ, 6, prefix='a') # or 4 or 5 or 11 ... a = A.basis() A.an_element() /// a[0] + 3*a[1] + 3*a[2] }}}

We now want to endow $A$ with its natural product structure, to get the desired group algebra. To define a multiplication, we should be in a category where multiplication makes sense, which is not yet the case:

{{{id=4| A.category() /// Category of modules with basis over Rational Field }}}

We can look at the available categories from the documentation in the reference manual: http://sagemath.com/doc/reference/categories.html

Or we can use introspection:

{{{id=5| sage.categories. # Look through the list of categories to pick one we want /// }}}

Once we have chosen an appropriate category (here :class:`AlgebrasWithBasis`), one can look at one example:

{{{id=6| E = AlgebrasWithBasis(QQ).example(); E /// An example of an algebra with basis: the free algebra on the generators ('a', 'b', 'c') over Rational Field }}} {{{id=7| e = E.an_element(); e /// B[word: ] + 2*B[word: a] + 3*B[word: b] }}}

and browse through its code:

{{{id=8| E?? /// }}}

This code is meant as a template from which to start implementing a new algebra. In particular it suggests that we need to implement the methods product_on_basis, one_basis, _repr_ and algebra_generators. Another way to get this list of methods is to ask the category (TODO: find a slicker idiom for this):

{{{id=9| from sage.misc.abstract_method import abstract_methods_of_class abstract_methods_of_class(AlgebrasWithBasis(QQ).element_class) /// {'required': [], 'optional': ['_add_', '_mul_']} }}} {{{id=10| abstract_methods_of_class(AlgebrasWithBasis(QQ).parent_class) /// {'required': ['__contains__'], 'optional': ['one_basis', 'product_on_basis']} }}}

Warning

the result above is not yet necessarily complete; many

required methods in the categories are not yet marked as :function:`abstract_methods`. We also recommend browsing the documentation of this category: :class:`AlgebrasWithBasis`.

Here is the obtained implementation of the group algebra:

{{{id=11| class MyCyclicGroupAlgebra(CombinatorialFreeModule): def __init__(self, R, n, **keywords): self._group = Zmod(n) CombinatorialFreeModule.__init__(self, R, self._group, category=AlgebrasWithBasis(R), **keywords) def product_on_basis(self, left, right): return self.monomial( left + right ) def one_basis(self): return self._group.zero() def algebra_generators(self): return Family( [self.monomial( self._group(1) ) ] ) def _repr_(self): return "Jason's group algebra of %s over %s"%(self._group, self.base_ring()) /// }}}

Some notes about this implementation:

Let us do some calculations:

{{{id=12| A = MyCyclicGroupAlgebra(QQ, 2, prefix='a') # or 4 or 5 or 11 ... a = A.basis(); f = A.an_element(); A, f /// (Jason's group algebra of the cyclic group Zmod(2) over Rational Field, a[0] + 3*a[1]) }}} {{{id=13| f * f /// 10*a[0] + 6*a[1] }}} {{{id=14| f. f.is_idempotent() /// False }}} {{{id=15| A.one() /// a[0] }}} {{{id=16| x = A.algebra_generators().first() # Typically x,y, ... = A.algebra_generators() [x^i for i in range(4)] /// [a[0], a[1], a[0], a[1]] }}} {{{id=17| g = 2*a[1]; (f + g)*f == f*f + g*f /// True }}}

This seems to work fine, but we would like to put more stress on our implementation to shake potential bugs out of it. To this end, we will use :class:`TestSuite`, a tool which will perform many routine tests on our algebra for us:

{{{id=18| TestSuite(A).run(verbose=True) /// running ._test_additive_associativity() . . . passrunning ._test_an_element() . . . passrunning ._test_associativity() . . . passrunning ._test_category() . . . passrunning ._test_elements() . . . Running the test suite of self.an_element() running ._test_category() . . . pass running ._test_not_implemented_methods() . . . pass running ._test_pickling() . . . pass passrunning ._test_not_implemented_methods() . . . passrunning ._test_one() . . . passrunning ._test_pickling() . . . passrunning ._test_prod() . . . passrunning ._test_some_elements() . . . passrunning ._test_zero() . . . pass }}}

For more information on categories, see:

{{{id=19| sage.categories.primer? /// }}}

Review

We wanted to create an algebra, so we:

# Created the underlying vector space using :class:`CombinatorialFreeModule`
# Looked at ``sage.categories.<tab>`` to find an appropriate category
# Loaded an example of that category to see what methods we needed to write
# Added the category information and other necessary methods to our class
# Ran :class:`TestSuite` to catch potential discrepancies

Exercises

  1. Make a tiny modification to product_on_basis in "MyCyclicGroupAlgebra" to implement the dual of the group algebra of the cyclic group instead of its group algebra (product given by $b_fb_g=delta_{f,g}bf$).

    Run the :class:`TestSuite` tests (you may ignore the "pickling" errors). What do you notice?

    Fix the implementation of one and check that the tests now pass.

    Add the hopf algebra structure. Hint: look at the example:

    {{{id=20| C = HopfAlgebrasWithBasis(QQ).example() /// }}}
  2. Given a set $S$, say:

    {{{id=21| S = Set([1,2,3,4,5]) /// }}}

    and a base ring, say:

    {{{id=22| R = QQ /// }}}

    implement an $R$-algebra:

    {{{id=23| F = SubsetAlgebraOnFundamentalBasis(S, R) # todo: not implemented /// }}}

    whose basis (b_s)_{s\subset S} is indexed by the subsets of S:

    {{{id=24| Subsets(S) /// Subsets of {1, 2, 3, 4, 5} }}}

    and where the product is defined by $b_s b_t = b_{scup t}$.

Morphisms

To better understand relationships between algebraic spaces, one wants to consider morphisms between them:

{{{id=25| A.module_morphism? A = MyCyclicGroupAlgebra(QQ, 2, prefix='a') B = MyCyclicGroupAlgebra(QQ, 6, prefix='b') A, B /// (Jason's group algebra of the cyclic group Zmod(2) over Rational Field, Jason's group algebra of the cyclic group Zmod(6) over Rational Field) }}} {{{id=26| def func_on_basis(g): r""" This function is the 'brain' of a (linear) morphism from A --> B. The input is the index of basis element of the domain (A). The output is an element of the codomain (B). """ if g==1: return B.monomial(Zmod(6)(3)) else: return B.one() /// }}}

We can now define a morphism which extends this function to $A$ by linearity:

{{{id=27| phi = A.module_morphism(func_on_basis, codomain=B) f = A.an_element() f /// a[0] + 3*a[1] }}} {{{id=28| phi(f) /// b[0] + 3*b[3] }}}

Exercise

Define a new free module In with basis indexed by the subsets of $S$, and a morphism phi from In to F defined by

.. math:: \phi(In_s) = \sum_{ t\subset s}  F_t



Diagonal and Triangular Morphisms

We now illustrate how to specify that a given morphism is diagonal or triangular with respect to some order on the basis which makes it invertible. Currently this feature requires the domain and codomain to have the same index set (in progress ...).

{{{id=29| X = CombinatorialFreeModule(QQ, Partitions(), prefix='x'); x = X.basis(); Y = CombinatorialFreeModule(QQ, Partitions(), prefix='y'); y = Y.basis(); /// }}}

A diagonal module morphism takes as argument a function whose input is the index of a basis element of the domain, and whose output is the coefficient of the corresponding basis element of the codomain:

{{{id=30| def diag_func(p): if len(p)==0: return 1 else: return p[0] diag_func(Partition([3,2,1])) /// 3 }}} {{{id=31| X_to_Y = X.module_morphism(diagonal=diag_func, codomain=Y) f = X.an_element(); f /// x[[]] + 2*x[[1]] + 3*x[[2]] }}} {{{id=32| X_to_Y(f) /// y[[]] + 2*y[[1]] + 6*y[[2]] }}}

Python fun-fact: ~ is the inversion operator (but be careful with int's!):

{{{id=33| ~2 /// 1/2 }}} {{{id=34| ~(int(2)) /// -3 }}}

Diagonal module morphisms are invertible:

{{{id=35| Y_to_X = ~X_to_Y f = y[Partition([3])] - 2*y[Partition([2,1])] f /// -2*y[[2, 1]] + y[[3]] }}} {{{id=36| Y_to_X(f) /// -x[[2, 1]] + 1/3*x[[3]] }}} {{{id=37| X_to_Y(Y_to_X(f)) /// -2*y[[2, 1]] + y[[3]] }}}

For triangular morphisms, just like ordinary morphisms, we need a function which accepts as input the index of a basis element of the domain and returns an element of the codomain. We think of this function as representing the columns of the matrix of the linear transformation:

{{{id=38| def triang_on_basis(p): return Y.sum_of_monomials(mu for mu in Partitions(sum(p)) if mu >= p) triang_on_basis([3,2]) /// y[[3, 2]] + y[[4, 1]] + y[[5]] }}} {{{id=39| X_to_Y = X.module_morphism(triang_on_basis, triangular='lower', unitriangular=True, codomain=Y) f = x[Partition([1,1,1])] + 2*x[Partition([3,2])]; f /// x[[1, 1, 1]] + 2*x[[3, 2]] }}} {{{id=40| X_to_Y(f) /// y[[1, 1, 1]] + y[[2, 1]] + y[[3]] + 2*y[[3, 2]] + 2*y[[4, 1]] + 2*y[[5]] }}}

Triangular module_morphisms are also invertible, even if X and Y are both infinite-dimensional:

{{{id=41| Y_to_X = ~X_to_Y f /// x[[1, 1, 1]] + 2*x[[3, 2]] }}} {{{id=42| Y_to_X(X_to_Y(f)) /// x[[1, 1, 1]] + 2*x[[3, 2]] }}}

For details, see :meth:`ModulesWithBasis.ParentMethods.module_morphism` (and also :class:`sage.categories.modules_with_basis.TriangularModuleMorphism`):

{{{id=43| A.module_morphism? /// }}}

Exercise

Redefine the morphism phi from the previous exercise to specify that it is triangular w.r.t. inclusion of subsets, and inverse this morphism. You may want to use the following comparison function as cmp argument to modules_morphism:

{{{id=44| def subset_cmp(s, t): """ A comparison function on sets which gives a linear extension of the inclusion order. INPUT: - ``x``, ``y`` -- sets EXAMPLES:: sage: sorted(Subsets([1,2,3]), subset_cmp) [{}, {1}, {2}, {3}, {1, 2}, {1, 3}, {2, 3}, {1, 2, 3}] """ s = cmp(len(x), len(y)) if s != 0: return s return cmp(list(x), list(y)) /// }}}

Coercions

Once we have defined a morphism from $X to Y$, we can register it as a coercion. This will allow Sage to apply the morphism automatically whenever we combine elements of $X$ and $Y$ together. See http://sagemath.com/doc/reference/coercion.html for more information. As a training step, let us first define a morphism $X$ to $h$, and register it as a coercion:

{{{id=45| def triang_on_basis(p): return h.sum_of_monomials(mu for mu in Partitions(sum(p)) if mu >= p) triang_on_basis([3,2]) /// h[3, 2] + h[4, 1] + h[5] }}} {{{id=46| X_to_h = X.module_morphism(triang_on_basis, triangular='lower', unitriangular=True, codomain=h) X_to_h. X_to_h.register_as_coercion() /// }}}

Now we can convert elements from $X$ to $h$, but also do mixed arithmetic with them:

{{{id=47| h(x[Partition([3,2])]) /// h[3, 2] + h[4, 1] + h[5] }}} {{{id=48| h([2,2,1]) + x[Partition([2,2,1])] /// 2*h[2, 2, 1] + h[3, 1, 1] + h[3, 2] + h[4, 1] + h[5] }}}

Exercise

Use the inverse of phi to implement the inverse coercion from F to In. Reimplement In as an algebra, with a product method making it use phi and its inverse.

Application: new basis and quotients of symmetric functions

In the sequel, we show how to combine everything we have seen to implement a new basis of the algebra of symmetric functions:

{{{id=49| SF = SymmetricFunctions(QQ); # A GradedHopfAlgebraWithBasis h = SF.homogeneous() # A particular basis, indexed by partitions (with some additional magic) /// }}}

$h$ is a graded algebra whose basis is indexed by partitions. Namely h([i]) represents the sum of all monomials of degree $i$:

sage: h([2]).expand(4) x0^2 + x0*x1 + x1^2 + x0*x2 + x1*x2 + x2^2 + x0*x3 + x1*x3 + x2*x3 + x3^2

and h(\mu) = prod( h(p) for p in mu ):

{{{id=50| h([3,2,2,1]) == h([3]) * h([2]) * h([2]) * h([1]) /// True }}} {{{id=51| class MySFBasis(CombinatorialFreeModule): r""" Note: We would typically use SymmetricFunctionAlgebra_generic for this. This is as an example only. """ def __init__(self, R, *args, **kwargs): """ TODO: Informative doc-string and examples """ CombinatorialFreeModule.__init__(self, R, Partitions(), category=AlgebrasWithBasis(R), *args, **kwargs) self._h = SymmetricFunctions(R).homogeneous() self._to_h = self.module_morphism( self._to_h_on_basis, triangular='lower', unitriangular=True, codomain=self._h) self._from_h = ~(self._to_h) self._to_h.register_as_coercion() self._from_h.register_as_coercion() def _to_h_on_basis(self, la): return self._h.sum_of_monomials(mu for mu in Partitions(sum(la)) if mu >= la) def product(self, left, right): return self( self._h(left) * self._h(right) ) def _repr_(self): return "Jason's basis for symmetric functions over %s"%self.base_ring() @cached_method def one_basis(self): r""" Returns the index of the basis element which is equal to '1'.""" return Partition([]) X = MySFBasis(QQ, prefix='x'); x = X.basis(); h = SymmetricFunctions(QQ).homogeneous() f = X(h([2,1,1]) - 2*h([2,2])) # Note the capital X f h(f) /// x[[2, 1, 1]] - 3*x[[2, 2]] + 2*x[[3, 1]]h[2, 1, 1] - 2*h[2, 2] }}} {{{id=52| f*f*f /// x[[2, 2, 2, 1, 1, 1, 1, 1, 1]] - 7*x[[2, 2, 2, 2, 1, 1, 1, 1]] + 18*x[[2, 2, 2, 2, 2, 1, 1]] - 20*x[[2, 2, 2, 2, 2, 2]] + 8*x[[3, 1, 1, 1, 1, 1, 1, 1, 1, 1]] }}} {{{id=53| h(f*f) /// h[2, 2, 1, 1, 1, 1] - 4*h[2, 2, 2, 1, 1] + 4*h[2, 2, 2, 2] }}}

As a final example, we implement a quotient of the algebra of symmetric functions:

{{{id=54| class MySFQuotient(CombinatorialFreeModule): r""" The quotient of the ring of symmetric functions by the ideal generated by those monomial symmetric functions whose part is larger than some fixed number ``k``. """ def __init__(self, R, k, prefix=None, *args, **kwargs): # Note: Setting self._prefix is equivalent to using the prefix keyword # in CombinatorialFreeModule.__init__ if prefix is not None: self._prefix = prefix else: self._prefix = 'mm' CombinatorialFreeModule.__init__(self, R, Partitions(NonNegativeIntegers(), max_part=k), category = GradedHopfAlgebrasWithBasis(R), *args, **kwargs) self._k = k self._m = SymmetricFunctions(R).monomial() self.lift = self.module_morphism(self._m.monomial) self.retract = self._m.module_morphism(self._retract_on_basis, codomain=self) self.lift.register_as_coercion() self.retract.register_as_coercion() def _retract_on_basis(self, mu): r""" Takes the index of a basis element of a monomial symmetric function, and returns the projection of that element to the quotient.""" if len(mu) > 0 and mu[0] > self._k: return self.zero() return self.monomial(mu) @cached_method def one_basis(self): return Partition([]) def product(self, left, right): return self( self._m(left) * self._m(right) ) MM = MySFQuotient(QQ, 3) mm = MM.basis() m = SymmetricFunctions(QQ).monomial() P = Partition f = mm[P([3,2,1])] + 2*mm[P([3,3])] f /// mm[[3, 2, 1]] + 2*mm[[3, 3]] }}} {{{id=55| m(f) /// m[3, 2, 1] + 2*m[3, 3] }}} {{{id=56| (m(f))^2 /// 8*m[3, 3, 2, 2, 1, 1] + 12*m[3, 3, 2, 2, 2] + 24*m[3, 3, 3, 2, 1] + 48*m[3, 3, 3, 3] + 4*m[4, 3, 2, 2, 1] + 4*m[4, 3, 3, 1, 1] + 14*m[4, 3, 3, 2] + 4*m[4, 4, 2, 2] + 4*m[4, 4, 3, 1] + 6*m[4, 4, 4] + 4*m[5, 3, 2, 1, 1] + 4*m[5, 3, 2, 2] + 12*m[5, 3, 3, 1] + 2*m[5, 4, 2, 1] + 6*m[5, 4, 3] + 4*m[5, 5, 1, 1] + 2*m[5, 5, 2] + 4*m[6, 2, 2, 1, 1] + 6*m[6, 2, 2, 2] + 6*m[6, 3, 2, 1] + 10*m[6, 3, 3] + 2*m[6, 4, 1, 1] + 5*m[6, 4, 2] + 4*m[6, 5, 1] + 4*m[6, 6] }}} {{{id=57| f^2 /// 8*mm[[3, 3, 2, 2, 1, 1]] + 12*mm[[3, 3, 2, 2, 2]] + 24*mm[[3, 3, 3, 2, 1]] + 48*mm[[3, 3, 3, 3]] }}} {{{id=58| (m(f))^2 - m(f^2) /// 4*m[4, 3, 2, 2, 1] + 4*m[4, 3, 3, 1, 1] + 14*m[4, 3, 3, 2] + 4*m[4, 4, 2, 2] + 4*m[4, 4, 3, 1] + 6*m[4, 4, 4] + 4*m[5, 3, 2, 1, 1] + 4*m[5, 3, 2, 2] + 12*m[5, 3, 3, 1] + 2*m[5, 4, 2, 1] + 6*m[5, 4, 3] + 4*m[5, 5, 1, 1] + 2*m[5, 5, 2] + 4*m[6, 2, 2, 1, 1] + 6*m[6, 2, 2, 2] + 6*m[6, 3, 2, 1] + 10*m[6, 3, 3] + 2*m[6, 4, 1, 1] + 5*m[6, 4, 2] + 4*m[6, 5, 1] + 4*m[6, 6] }}} {{{id=59| MM( (m(f))^2 - m(f^2) ) /// 0 }}}