# -*- coding: utf-8 -*-
r"""
Base class for the Geometric Representations of Coxeter Groups
A geometric representations of a Coxeter group is a morphism from the Coxeter Group to
a reflection group stabilizing a (sometimes non-canonical) bilinear form obtained from
the group. This class allows to do computations related to the geometric
representation. It is possible to compute roots, weights, images of the
elements, it is also possible to visualize several related geometric objects
such as the Tits cone, the isotropic cone, the limit roots, and more.
EXAMPLES::
sage: from brocoli import *
sage: M1 = CoxeterMatrix([[1,4,4],[4,1,4],[4,4,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1);GR1
Geometric representation of a Coxeter group of rank 3 with Coxeter matrix
[1 4 4]
[4 1 4]
[4 4 1]
If the parameter ``exact`` is set to ``False``, then the representation is done in
``RDF``::
sage: GR1rdf = GeometricRepresentationCoxeterGroup(M1,exact=False)
sage: GR1rdf.base_ring()
Real Double Field
sage: GR1.base_ring()
Universal Cyclotomic Field
It is possible to specify a the symbols (or alphabet) for the generators::
sage: M2 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,oo,oo],[oo,oo,1,oo],[oo,oo,oo,1]])
sage: GR2 = GeometricRepresentationCoxeterGroup(M2, generators=['a','b','c','d'])
sage: GR2.generators()
('a', 'b', 'c', 'd')
The geometric representation uses a certain bilinear form::
sage: GR1.bilinear_form()
[ 1 -1/2*E(8) + 1/2*E(8)^3 -1/2*E(8) + 1/2*E(8)^3]
[-1/2*E(8) + 1/2*E(8)^3 1 -1/2*E(8) + 1/2*E(8)^3]
[-1/2*E(8) + 1/2*E(8)^3 -1/2*E(8) + 1/2*E(8)^3 1]
sage: GR2.bilinear_form()
[ 1 -1 -1 -1]
[-1 1 -1 -1]
[-1 -1 1 -1]
[-1 -1 -1 1]
We can ask for the signature of the bilinear form::
sage: GR1.signature()
(2, 1, 0)
sage: GR2.signature()
(3, 1, 0)
so both of them are Lorentzian; they act on Lorentz space::
sage: GR2.is_lorentzian()
True
We can go from generators to their matrices::
sage: GR2.generator_to_reflection('a')
[-1 2 2 2]
[ 0 1 0 0]
[ 0 0 1 0]
[ 0 0 0 1]
sage: GR2.simple_reflections()[0]
[-1 2 2 2]
[ 0 1 0 0]
[ 0 0 1 0]
[ 0 0 0 1]
One can compute some elements of length two (the first two)::
sage: GR2.elements(2)[:2]
[
[ 1 0 0 0] [ 3 6 6 -2]
[ 2 -1 2 2] [ 0 1 0 0]
[ 0 0 1 0] [ 0 0 1 0]
[ 6 -2 6 3], [ 2 2 2 -1]
]
We can then ask for a reduced expression for an element that was obtained::
sage: GR2.matrix_to_word(GR2.elements(2)[0])
word: bd
sage: GR2.matrix_to_word(GR2.elements(2)[1])
word: da
We can compute some roots::
sage: GR1.roots(1)
{(1, E(8) - E(8)^3, 0), (1, 0, E(8) - E(8)^3), (E(8) - E(8)^3, 1, 0), (E(8) - E(8)^3, 0, 1), (0, E(8) - E(8)^3, 1), (0, 1, E(8) - E(8)^3)}
sage: GR2.roots(1)[:5]
[(0, 2, 0, 1), (0, 0, 2, 1), (0, 2, 1, 0), (0, 1, 0, 2), (2, 0, 0, 1)]
and some weights::
sage: GR1.fundamental_weights()
((1 - E(8) + E(8)^3, -1, -1),
(-1, 1 - E(8) + E(8)^3, -1),
(-1, -1, 1 - E(8) + E(8)^3))
sage: GR2.fundamental_weights()
((1/4, -1/4, -1/4, -1/4),
(-1/4, 1/4, -1/4, -1/4),
(-1/4, -1/4, 1/4, -1/4),
(-1/4, -1/4, -1/4, 1/4))
Finally, we can visualize the isotropic cone::
sage: GR1.visualize_isotropic_cone()
Graphics object consisting of 4 graphics primitives
sage: GR2.visualize_isotropic_cone()
Graphics3d Object
But that's a bit boring, so we can add roots::
sage: GR1.visualize_roots([0,1,2])
Graphics object consisting of 22 graphics primitives
sage: GR2.visualize_roots([0,1,2])
Graphics3d Object
or weights::
sage: GR1.visualize_weights([0,1,2])
Graphics object consisting of 19 graphics primitives
sage: GR2.visualize_weights([0,1,2])
Graphics3d Object
or limit roots::
sage: GR1.visualize_limit_roots([3,4,5])
Graphics object consisting of 51 graphics primitives
sage: GR2.visualize_limit_roots([2,3,4])
Graphics3d Object
we can add the isotropic cone::
sage: GR1.visualize_limit_roots([3,4,5],isotropic=True)
Graphics object consisting of 55 graphics primitives
sage: GR2.visualize_limit_roots([2,3,4],isotropic=True)
Graphics3d Object
and finally the Tits cone::
sage: GR1.visualize_tits_cone(2)
Graphics object consisting of 13 graphics primitives
sage: GR2.visualize_tits_cone(3)
Graphics3d Object
REFERENCES:
- Chapter 5 of Reflection Groups and Coxeter Groups, J. Humphreys, (1990)
- Discrete linear groups generated by reflections (Russian translated to English),
E.B. Vinberg, (1971)
- C. Hohlweg, J.-P. Labbé, and V. Ripoll, Asymptotical behaviour of roots of
infinite Coxeter groups, Canad. J. Math., **66**, (2014), no. 2, 323–353.
- M. Dyer, C. Hohlweg, V. Ripoll, Imaginary cones and limit roots of
infinite Coxeter groups, Math. Z. **284** (2016), no. 3-4, 715–780.
- C. Hohlweg, J.-P. Préaux, V. Ripoll, On the Limit Set of Root Systems of
Coxeter Groups acting on Lorentzian spaces, arxiv:1305.0052 (July 2013),
24 p.
- H. Chen and J.-P. Labbé, Lorentzian Coxeter systems and Boyd–Maxwell ball
packings, Geom. Dedic., **174**, (2015), no. 1, 43–73.
- H. Chen and J.-P. Labbé, Limit directions for Lorentzian Coxeter systems,
Groups Geom. Dyn. **11** (2017), no. 2, 469–498.
- C. Hohlweg and J.-P. Labbé, On inversion sets and the weak order in
Coxeter groups, European J. Combin. **55** (2016), 1--19.
AUTHORS:
- Jean-Philippe Labbé (2011-): Initial version
- Vivien Ripoll (2011-2013): Added more options and functionalities
"""
##############################################################################
# Copyright (C) 2011 Jean-Philippe Labbe <labbe at math.fu-berlin.de>
# 2011 Vivien Ripoll < vivien.ripoll at univie.ac.at>
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 2 of the License, or
# (at your option) any later version.
# http://www.gnu.org/licenses/
##############################################################################
from sage.arith.srange import xsrange
from sage.misc.cachefunc import cached_method
from sage.misc.flatten import flatten
from sage.misc.functional import n
from sage.misc.latex import LatexExpr
from sage.rings.rational_field import QQ
from sage.rings.qqbar import QQbar
from sage.rings.real_double import RDF
from sage.rings.real_mpfr import RR
from sage.rings.complex_double import CDF
from sage.rings.polynomial.polynomial_ring_constructor import PolynomialRing
from sage.rings.universal_cyclotomic_field import UniversalCyclotomicField
from sage.rings.infinity import Infinity as oo
from sage.calculus.var import var
from sage.symbolic.all import i as I
from sage.symbolic.relation import solve
from sage.functions.other import sqrt
from sage.functions.generalized import sign
from sage.functions.trig import arctan
from math import pi
from sage.matrix.constructor import matrix
from sage.modules.free_module_element import vector
from sage.modules.free_module import VectorSpace
from sage.geometry.polyhedron.constructor import Polyhedron
from sage.geometry.polyhedron.library import polytopes
from sage.categories.rings import Rings
from sage.combinat.words.word import Word
from sage.sets.set import Set
from sage.plot.point import point
from sage.plot.contour_plot import implicit_plot
from sage.plot.graphics import Graphics
from sage.plot.circle import circle
from sage.plot.arrow import arrow
from sage.plot.line import line
from sage.plot.plot3d.implicit_plot3d import implicit_plot3d
[docs]def reflection_image(reflection_vector, element, bilin_form_matrix):
"""
Returns the image of the vector ``element`` by the reflection
through the vector `reflection_vector` with respect to the bilinear form described
by the matrix ``bilin_form_matrix``.
INPUT:
- ``reflection_vector`` -- vector; a non-isotropic (pseudonorm 1)-vector to reflect with
- ``element`` -- vector; a vector that will be reflected
- ``bilin_form_matrix`` -- matrix; the matrix giving the bilinear form
OUTPUT:
The reflected vector
EXAMPLES::
sage: from brocoli import *
sage: bf = matrix([[1,-2],[-2,1]])
sage: reflection_image(vector([1,0]),vector([10,1]),bf)
(-6, 1)
sage: bf = matrix([[1,-1],[-1,1]])
sage: reflection_image(vector([1,0]),vector([1,1]),bf)
(1, 1)
"""
double_projection = (2 * reflection_vector * bilin_form_matrix * (element.column()))[0]
# pseudonorm = (reflection_vector*bilin_form_matrix*(reflection_vector.column()))[0]
# return element-(double_projection/pseudonorm)*reflection_vector
return element - double_projection * reflection_vector
[docs]def affinely_project_vector(vect, projection_space, affine_basis):
"""
This function returns the affine normalization of the vector ``vect`` in the
affine hyperplane orthogonal to the vector (1,...,1) according to the
affine basis ``affine_basis`` of the projection space ``projection_space``.
First, the vector is normalize to be a vector with sum of the coordinate
equals to 1. If the sum of the coordinates is equal to 0, it returns an
error. Then, the coefficients are used to give a linear combination of the
``affine_basis`` in the affine space `projection_space`.
For example, a 3-dimensional vector with have 3 coefficients which will be
used to obtain a linear combination of 3 vectors, say in a 2-dimensional
space.
INPUT:
- ``vect`` -- vector; the vector that is to project
- ``projection_space`` -- vector space; the affine space where to project
- ``affine_basis`` -- list of vectors; the image of the basis of the vector space of ``vect`` under a
projection to the affine space ``projection_space``
OUTPUT:
A vector of the affine space ``projection_space``
EXAMPLES::
sage: from brocoli import *
sage: PS = VectorSpace(RDF, 2)
sage: ab = regular_simplex_vertices(2)
sage: v = vector([1,2,3])
sage: affinely_project_vector(v,PS,ab)
(1.1666666666666665, 0.8660254037844386)
sage: UCF.<E> = UniversalCyclotomicField()
sage: u = vector(UCF, [2 + E(8) - E(8)^3, 1, E(8) - E(8)^3])
sage: affinely_project_vector(u,PS,ab)
(0.585786437626905, 0.4202659980740251)
sage: v = vector([1,-4,3])
sage: affinely_project_vector(v,PS,ab)
Traceback (most recent call last):
...
ValueError: The vector (1, -4, 3) does not have an affine image in the affine basis
"""
dim = len(affine_basis)
height = RDF(sum(vect))
if height == 0:
raise ValueError("The vector {} does not have an affine image in the affine basis".format(vect.__repr__()))
Nvect = vector([RDF(i)/height for i in vect])
image = projection_space(sum(Nvect[j]*affine_basis[j] for j in range(dim)))
return image
[docs]def regular_simplex_vertices(dim):
r"""
This function returns the floating point approximation of vertices of a regular
simplex of dimension ``dim`` having ``dim+1`` vertices
INPUT:
- ``dim`` : integer; the dimension of the desired simplex
OUTPUT:
A list of vectors, the vertices of a regular simplex
EXAMPLES::
sage: from brocoli import *
sage: regular_simplex_vertices(1)
[(0.0, 0.0), (1.0, 0.0)]
sage: regular_simplex_vertices(2)
[(0.0, 0.0), (2.0, 0.0), (1.0, 1.7320508075688772)]
sage: regular_simplex_vertices(3)
[(0.0, 0.0, 0.0),
(2.0, 0.0, 0.0),
(1.0, 1.7320508075688772, 0.0),
(1.0, 0.5773502691896257, 1.6329931618554518)]
sage: regular_simplex_vertices(4)
Traceback (most recent call last):
...
NotImplementedError: dimension >=4
"""
if dim == 1:
vs = VectorSpace(RDF, 2)
return [vs([0, 0]), vs([1, 0])]
elif dim == 2:
vs = VectorSpace(RDF, dim)
return [vs([0, 0]), vs([2, 0]), vs([1, sqrt(3)])]
elif dim == 3:
vs = VectorSpace(RDF, dim)
return [vs([0, 0, 0]), vs([2, 0, 0]), vs([1, sqrt(3), 0]),
vs([1, 1/sqrt(3), (2*sqrt(2))/sqrt(3)])]
else:
raise NotImplementedError("dimension >=4")
[docs]def plot_simplex(size, color=(0, 1, 0)):
"""
This function returns a graphical element representing the regular simplex
in dimension 2 or 3
INPUT:
- ``size`` -- integer; the number of vertices of the simplex
- ``color`` -- a rgb vector; the color of the simplex
OUTPUT:
A graphical object
EXAMPLES::
sage: from brocoli import *
sage: plot_simplex(2)
Graphics object consisting of 1 graphics primitive
sage: plot_simplex(3)
Graphics object consisting of 3 graphics primitives
sage: plot_simplex(5)
Traceback (most recent call last):
...
ValueError: The size of the simplex should be 2, 3 or 4
"""
if size == 4:
simplex_image = Graphics()
projection_space = VectorSpace(RDF, size-1)
affine_basis = regular_simplex_vertices(size-1)
canonical_basis = [vector([1, 0, 0, 0]),
vector([0, 1, 0, 0]),
vector([0, 0, 1, 0]),
vector([0, 0, 0, 1])]
simplex_vertices = [affinely_project_vector(pt, projection_space, affine_basis) for pt in canonical_basis]
for i in range(4):
for j in range(4):
if i < j: # Adding the edges of the simplex
simplex_image += line([simplex_vertices[i], simplex_vertices[j]], color=color, thickness=2)
elif size == 3: # Adding the edges of the simplex
simplex_image = line([(0, 0), (1, sqrt(3))], color=color) + \
line([(0, 0), (2, 0)], color=color) + \
line([(2, 0), (1, sqrt(3))], color=color)
elif size == 2:
from sage.plot.point import points
simplex_image = points([(1, 0), (0, 1)], color=color)
else:
raise ValueError("The size of the simplex should be 2, 3 or 4")
return simplex_image
[docs]class GeometricRepresentationCoxeterGroup():
r"""
Base class for geometric representations of Coxeter groups.
INPUT:
- ``coxeter_matrix`` -- a CoxeterMatrix object form Sage
- ``generators`` -- an alphabet for the generators (default: ``None``)
- ``exact`` -- boolean keyword argument (default: ``True``); whether to do the
computation in an exact ring.
EXAMPLES:
To create a geometric representation, we start with Coxeter matrices::
sage: from brocoli import *
sage: QF.<a> = QuadraticField(2)
sage: M1 = CoxeterMatrix([[1,4,4],[4,1,4],[4,4,1]])
sage: M2 = CoxeterMatrix([[1,4,oo],[4,1,4],[oo,4,1]])
sage: M3 = CoxeterMatrix([[1,4,-2],[4,1,4],[-2,4,1]])
sage: M4 = CoxeterMatrix([[1,4,-1.5],[4,1,4],[-1.5,4,1]])
sage: M5 = CoxeterMatrix([[1,4,-3/2],[4,1,4],[-3/2,4,1]])
sage: M6 = CoxeterMatrix([[1,4,-3/2],[4,1,-a],[-3/2,-a,1]])
sage: M7 = CoxeterMatrix([[1,4,-1.5],[4,1,-a],[-1.5,-a,1]])
depending on the Coxeter matrix, different base rings are chosen for the
representation. When some labels are integers greater or equal to 4, then
it uses the Universal Cyclotomic Field. If some other not rational labels
are given, then it tries to find the appropriate ring::
sage: for m in [M1,M2,M3,M4,M5,M6,M7]:
....: print(GeometricRepresentationCoxeterGroup(m).base_ring())
....:
Universal Cyclotomic Field
Universal Cyclotomic Field
Universal Cyclotomic Field
Real Double Field
Universal Cyclotomic Field
Algebraic Field
Real Double Field
There is a slightly different behavior if there are no labels that require
the Univeral Cyclotomic Field::
sage: M1p = CoxeterMatrix([[1,3,3],[3,1,3],[3,3,1]])
sage: M2p = CoxeterMatrix([[1,3,oo],[3,1,3],[oo,3,1]])
sage: M3p = CoxeterMatrix([[1,3,-2],[3,1,3],[-2,3,1]])
sage: M4p = CoxeterMatrix([[1,3,-1.5],[3,1,3],[-1.5,3,1]])
sage: M5p = CoxeterMatrix([[1,3,-3/2],[3,1,3],[-3/2,3,1]])
sage: M6p = CoxeterMatrix([[1,3,-3/2],[3,1,-a],[-3/2,-a,1]])
sage: M7p = CoxeterMatrix([[1,3,-1.5],[3,1,-a],[-1.5,-a,1]])
sage: for m in [M1p,M2p,M3p,M4p,M5p,M6p,M7p]:
....: print(GeometricRepresentationCoxeterGroup(m).base_ring())
....:
Rational Field
Rational Field
Rational Field
Real Double Field
Rational Field
Number Field in a with defining polynomial x^2 - 2
Real Double Field
If the parameter ``exact`` is set to ``False``, then computations are done
in ``RDF``::
sage: GR1rdf = GeometricRepresentationCoxeterGroup(M1,exact=False);GR1rdf.base_ring()
Real Double Field
It is possible to specify a the symbols (or alphabet) for the generators::
sage: GR1 = GeometricRepresentationCoxeterGroup(M1, generators=['a','b','c'])
sage: GR1.generators()
('a', 'b', 'c')
TESTS::
sage: S = [GeometricRepresentationCoxeterGroup(CM) for CM in CoxeterMatrix.samples()]
"""
def __init__(self, coxeter_matrix, generators=None, exact=True):
"""
Geometric representation of a Coxeter group constructor.
TESTS::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,oo,4],[oo,1,oo],[4,oo,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: TestSuite(GR).run(skip='_test_pickling')
"""
# Setting up the general data structures
self._rank = coxeter_matrix.rank()
self._coxeter_matrix_obj = coxeter_matrix
self._coxeter_matrix = self._coxeter_matrix_obj._matrix
if generators is None:
self._generators = tuple(range(1, self._rank+1))
elif len(generators) != self._rank:
raise ValueError("The number of generators is not equal to the rank")
else:
self._generators = tuple(generators)
# Initializing the set of computed limit roots
self._computed_limit_roots = Set()
# Recording the labels of the Coxeter graph
labels = list(set(self._coxeter_matrix.list()))
labels.sort()
labels.remove(1)
self._edge_labels = tuple(labels)
# Check if need the usage of cyclotomic field
need_cyclotomic = any(lab not in [oo, 2, 3] and lab > 0 for lab in labels)
# Getting the labels that are potentially irrational
potential_weird_infinite_values = [lab for lab in labels if lab < 0 and not QQ.has_coerce_map_from(lab.parent())]
# Setting up the base ring
if not exact:
# Computations done in RDF
base_ring = RDF
elif any(x.parent() in [RDF, RR] for x in potential_weird_infinite_values):
# Force RDF when found in input
base_ring = RDF
elif not need_cyclotomic:
# No labels >= 4
if not potential_weird_infinite_values:
# No entries leading to irrational values
base_ring = QQ
else:
# We analyze the given labels to find the correct base_ring
from sage.structure.sequence import Sequence
from sage.structure.coerce import py_scalar_parent
the_parent = Sequence(potential_weird_infinite_values).universe()
if isinstance(the_parent, type):
base_ring = py_scalar_parent(the_parent)
else:
base_ring = the_parent
if base_ring is RR:
base_ring = RDF
if base_ring not in Rings() or not RDF.has_coerce_map_from(base_ring):
raise ValueError("the input values could not lead to a common ring.")
elif not potential_weird_infinite_values:
# No entries leading to irrational values, so we use UniversalCyclotomicField
base_ring = UniversalCyclotomicField()
else:
# There are labels >= 4 and potential weird infinite entries
base_ring = QQbar # Find a better ring for more specific entries
self._base_ring = base_ring
# Setting up dictionaries
self._matrix_to_word = {}
self._generator_to_reflection = {}
# Setting up computational sets
self._computed_roots = Set()
self._computed_weights_exact = Set()
self._computed_space_weights_exact = Set()
self._fund_weights_matrix_exact = matrix()
self._computed_weights_rdf = Set()
self._computed_space_weights_rdf = Set()
self._fund_weights_matrix_rdf = matrix()
def __repr__(self):
r"""
Return a string representation of ``self``.
EXAMPLES::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,oo,oo],[oo,1,oo],[oo,oo,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M);GR
Geometric representation of a Coxeter group of rank 3 with Coxeter matrix
[ 1 -1 -1]
[-1 1 -1]
[-1 -1 1]
"""
msg = 'Geometric representation of a Coxeter group of rank {} with Coxeter matrix\n{}'
return msg.format(self._rank, self._coxeter_matrix)
[docs] def base_ring(self):
r"""
Return the base ring of the representation
OUTPUT:
A ring
EXAMPLES::
sage: from brocoli import *
sage: QF.<a> = QuadraticField(2)
sage: M1 = CoxeterMatrix([[1,oo,oo],[oo,1,oo],[oo,oo,1]])
sage: M2 = CoxeterMatrix([[1,-3/2,-3/2],[-3/2,1,-3/2],[-3/2,-3/2,1]])
sage: M3 = CoxeterMatrix([[1,-1.5,-1.5],[-1.5,1,-1.5],[-1.5,-1.5,1]])
sage: M4 = CoxeterMatrix([[1,3,4],[3,1,5],[4,5,1]])
sage: M5 = CoxeterMatrix([[1,3,-3/2],[3,1,-a],[-3/2,-a,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1);GR1.base_ring()
Rational Field
sage: GR2 = GeometricRepresentationCoxeterGroup(M2);GR2.base_ring()
Rational Field
sage: GR3 = GeometricRepresentationCoxeterGroup(M3);GR3.base_ring()
Real Double Field
sage: GR4 = GeometricRepresentationCoxeterGroup(M4);GR4.base_ring()
Universal Cyclotomic Field
sage: GR5 = GeometricRepresentationCoxeterGroup(M5);GR5.base_ring()
Number Field in a with defining polynomial x^2 - 2
"""
return self._base_ring
[docs] def coxeter_matrix(self):
r"""
Return the associated Coxeter matrix
OUTPUT:
A Coxeter Matrix object
EXAMPLES::
sage: from brocoli import *
sage: M1 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,oo,oo],[oo,oo,1,oo],[oo,oo,oo,1]])
sage: M2 = CoxeterMatrix([[1,3,3,3],[3,1,3,3],[3,3,1,3],[3,3,3,1]])
sage: M3 = CoxeterMatrix([[1,4,-2],[4,1,4],[-2,4,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: GR2 = GeometricRepresentationCoxeterGroup(M2)
sage: GR3 = GeometricRepresentationCoxeterGroup(M3)
sage: GR1.coxeter_matrix()
[ 1 -1 -1 -1]
[-1 1 -1 -1]
[-1 -1 1 -1]
[-1 -1 -1 1]
sage: GR2.coxeter_matrix()
[1 3 3 3]
[3 1 3 3]
[3 3 1 3]
[3 3 3 1]
sage: GR3.coxeter_matrix()
[ 1 4 -2]
[ 4 1 4]
[-2 4 1]
"""
return self._coxeter_matrix_obj
[docs] def edge_labels(self):
r"""
Return the edge labels of the associated the Coxeter graph
OUTPUT:
a tuple containing the labels
EXAMPLES::
sage: from brocoli import *
sage: M1 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,oo,oo],[oo,oo,1,oo],[oo,oo,oo,1]])
sage: M2 = CoxeterMatrix([[1,3,3,3],[3,1,3,-1.5],[3,3,1,-3/2],[3,-1.5,-3/2,1]])
sage: M3 = CoxeterMatrix([[1,4,-2],[4,1,3],[-2,3,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: GR2 = GeometricRepresentationCoxeterGroup(M2)
sage: GR3 = GeometricRepresentationCoxeterGroup(M3)
sage: GR1.edge_labels()
(-1,)
sage: GR2.edge_labels()
(-1.50000000000000, 3.00000000000000)
sage: GR3.edge_labels()
(-2, 3, 4)
"""
return self._edge_labels
[docs] def generators(self):
r"""
Return the generators of the associated Coxeter group as letters of an alphabet
OUTPUT:
A tuple containing the (letter)-generators of the Coxeter group
EXAMPLES::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,oo,oo,oo],[oo,1,oo,oo],[oo,oo,1,oo],[oo,oo,oo,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M)
sage: GR1.generators()
(1, 2, 3, 4)
sage: GR2 = GeometricRepresentationCoxeterGroup(M,generators=['a','b','c','d'])
sage: GR2.generators()
('a', 'b', 'c', 'd')
sage: GR3 = GeometricRepresentationCoxeterGroup(M,generators=['a','b','c'])
Traceback (most recent call last):
...
ValueError: The number of generators is not equal to the rank
"""
return self._generators
[docs] def rank(self):
r"""
Return the rank of the associated Coxeter group
OUTPUT:
An integer
EXAMPLES::
sage: from brocoli import *
sage: B2 = CoxeterMatrix([[1,4],[4,1]])
sage: GR_B2 = GeometricRepresentationCoxeterGroup(B2)
sage: GR_B2.rank()
2
sage: A3 = CoxeterMatrix([[1,3,2],[3,1,3],[2,3,1]])
sage: GR_A3 = GeometricRepresentationCoxeterGroup(A3)
sage: A3.rank()
3
sage: M = CoxeterMatrix([[1,oo,oo,oo],[oo,1,oo,oo],[oo,oo,1,oo],[oo,oo,oo,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: GR.rank()
4
"""
return self._rank
[docs] @cached_method
def identity_element(self):
r"""
Return the matrix corresponding to the identity element, i.e. the
identity matrix of the representation space.
OUTPUT:
an identity matrix
EXAMPLES::
sage: from brocoli import *
sage: M1 = CoxeterMatrix([[1,4,4],[4,1,4],[4,4,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: GR1.identity_element()
[1 0 0]
[0 1 0]
[0 0 1]
sage: GR1.identity_element().base_ring()
Universal Cyclotomic Field
sage: M2 = CoxeterMatrix([[1,4,4,-2],[4,1,4,-1],[4,4,1,-1.5],[-2,-1,-1.5,1]])
sage: GR2 = GeometricRepresentationCoxeterGroup(M2)
sage: GR2.identity_element()
[1.0 0.0 0.0 0.0]
[0.0 1.0 0.0 0.0]
[0.0 0.0 1.0 0.0]
[0.0 0.0 0.0 1.0]
sage: GR2.identity_element().base_ring()
Real Double Field
"""
from sage.matrix.special import identity_matrix
identity = identity_matrix(self._base_ring, self._rank)
identity.set_immutable()
self._matrix_to_word[identity] = Word([], self._generators)
return identity
[docs] @cached_method
def simple_roots(self):
r"""
Return the simple roots
OUTPUT:
A tuple containing the simple roots as vectors
EXAMPLES::
sage: from brocoli import *
sage: M1 = CoxeterMatrix([[1,4,4],[4,1,4],[4,4,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: GR1.simple_roots()
((1, 0, 0), (0, 1, 0), (0, 0, 1))
sage: M2 = CoxeterMatrix([[1,4,4,-2],[4,1,4,-1],[4,4,1,-1.5],[-2,-1,-1.5,1]])
sage: GR2 = GeometricRepresentationCoxeterGroup(M2)
sage: GR2.simple_roots()
((1.0, 0.0, 0.0, 0.0),
(0.0, 1.0, 0.0, 0.0),
(0.0, 0.0, 1.0, 0.0),
(0.0, 0.0, 0.0, 1.0))
"""
simple_roots = [vector(self._base_ring, [0]*i + [1] + [0]*(self._rank - i - 1)) for i in range(self._rank)]
for vect in simple_roots:
vect.set_immutable()
self._computed_roots = Set(simple_roots)
return tuple(simple_roots)
[docs] @cached_method
def simple_reflections(self):
r"""
Return the reflections associated to the simple roots
OUTPUT:
a tuple containing the simple reflections as matrices
EXAMPLES::
sage: from brocoli import *
sage: M1 = CoxeterMatrix([[1,4,4],[4,1,4],[4,4,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: GR1.simple_reflections()[1]
[ 1 0 0]
[E(8) - E(8)^3 -1 E(8) - E(8)^3]
[ 0 0 1]
sage: M2 = CoxeterMatrix([[1,4,4,-2],[4,1,4,-1],[4,4,1,-1.5],[-2,-1,-1.5,1]])
sage: GR2 = GeometricRepresentationCoxeterGroup(M2)
sage: GR2.simple_reflections()[3]
[ 1.0 0.0 0.0 0.0]
[ 0.0 1.0 0.0 0.0]
[ 0.0 0.0 1.0 0.0]
[ 4.0 2.0 3.0 -1.0]
"""
simple_reflections = []
for index in range(self._rank):
s_root = self.simple_roots()[index]
generator = self._generators[index]
columns = [reflection_image(s_root, other_root, self.bilinear_form()) for other_root in self.simple_roots()]
new_matrix = matrix(self._base_ring, columns).transpose()
simple_reflections += [new_matrix]
simple_reflections[-1].set_immutable()
self._matrix_to_word[simple_reflections[-1]] = Word([generator], self._generators)
self._generator_to_reflection[generator] = simple_reflections[-1]
return tuple(simple_reflections)
[docs] def generator_to_reflection(self, generator):
r"""
Return the matrix associated to ``generator``
INPUT:
- ``generator`` -- letter or number; a generator of the Coxeter group
OUTPUT:
a matrix representing the generator
EXAMPLES::
sage: from brocoli import *
sage: M1 = CoxeterMatrix([[1,4,4],[4,1,4],[4,4,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: GR1.generator_to_reflection(1)
[ -1 E(8) - E(8)^3 E(8) - E(8)^3]
[ 0 1 0]
[ 0 0 1]
sage: GR1.generator_to_reflection(2)
[ 1 0 0]
[E(8) - E(8)^3 -1 E(8) - E(8)^3]
[ 0 0 1]
sage: GR1.generator_to_reflection(3)
[ 1 0 0]
[ 0 1 0]
[E(8) - E(8)^3 E(8) - E(8)^3 -1]
Also works if we change the set of generators::
sage: M2 = CoxeterMatrix([[1,oo,oo],[oo,1,oo],[oo,oo,1]])
sage: GR2 = GeometricRepresentationCoxeterGroup(M2,generators=['a','b','c'])
sage: GR2.generator_to_reflection(1)
Traceback (most recent call last):
...
KeyError: 1
sage: GR2.generator_to_reflection('a')
[-1 2 2]
[ 0 1 0]
[ 0 0 1]
sage: GR2.generator_to_reflection('b')
[ 1 0 0]
[ 2 -1 2]
[ 0 0 1]
sage: GR2.generator_to_reflection('c')
[ 1 0 0]
[ 0 1 0]
[ 2 2 -1]
"""
self.simple_reflections() # to make sure the dictionary is filled
return self._generator_to_reflection[generator]
[docs] @cached_method
def fundamental_weights(self, exact=True):
r"""
Return the fundamental weights of the representation
INPUT:
- ``exact`` -- a boolean (default = ``True``); whether to give the bilinear in an exact base
ring or not.
OUTPUT:
A tuple of vectors
EXAMPLES::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,oo,oo,oo],[oo,1,oo,oo],[oo,oo,1,oo],[oo,oo,oo,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: GR.fundamental_weights()
((1/4, -1/4, -1/4, -1/4),
(-1/4, 1/4, -1/4, -1/4),
(-1/4, -1/4, 1/4, -1/4),
(-1/4, -1/4, -1/4, 1/4))
sage: GR.fundamental_weights(False)
((0.25, -0.25, -0.25, -0.25),
(-0.25, 0.25, -0.25, -0.25),
(-0.25, -0.25, 0.25, -0.25),
(-0.25, -0.25, -0.25, 0.25))
TESTS::
sage: N = CoxeterMatrix([[1,oo,oo,oo],[oo,1,oo,oo],[oo,oo,1,oo],[oo,oo,oo,1]])
sage: GRN = GeometricRepresentationCoxeterGroup(N)
sage: GRN._computed_weights_exact
{}
sage: GRN.fundamental_weights()
((1/4, -1/4, -1/4, -1/4),
(-1/4, 1/4, -1/4, -1/4),
(-1/4, -1/4, 1/4, -1/4),
(-1/4, -1/4, -1/4, 1/4))
sage: GRN._computed_weights_exact
{(-1/4, 1/4, -1/4, -1/4), (-1/4, -1/4, 1/4, -1/4), (1/4, -1/4, -1/4, -1/4), (-1/4, -1/4, -1/4, 1/4)}
sage: GRN._computed_weights_rdf
{}
sage: GRN.fundamental_weights(False)
((0.25, -0.25, -0.25, -0.25),
(-0.25, 0.25, -0.25, -0.25),
(-0.25, -0.25, 0.25, -0.25),
(-0.25, -0.25, -0.25, 0.25))
sage: GRN._computed_weights_rdf
{(0.25, -0.25, -0.25, -0.25), (-0.25, 0.25, -0.25, -0.25), (-0.25, -0.25, 0.25, -0.25), (-0.25, -0.25, -0.25, 0.25)}
"""
bf = self.bilinear_form(exact)
if bf.det() != 0:
fundamental_weights = bf.inverse().columns()
for fw in fundamental_weights:
fw.set_immutable()
if exact:
self._computed_weights_exact = Set(fundamental_weights)
else:
self._computed_weights_rdf = Set(fundamental_weights)
else:
fundamental_weights = []
if exact:
self._fund_weights_matrix_exact = matrix(fundamental_weights).transpose()
else:
self._fund_weights_matrix_rdf = matrix(fundamental_weights).transpose()
return tuple(fundamental_weights)
[docs] @cached_method
def fundamental_space_weights(self, exact=True):
r"""
Return the fundamental weights which are positive with respect to the bilinear form.
INPUT:
- ``exact`` -- a boolean (default = ``True``); whether to give the bilinear in an exact base
ring or not.
OUTPUT:
a list of vectors
EXAMPLES::
sage: from brocoli import *
sage: M1 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,oo,oo],[oo,oo,1,oo],[oo,oo,oo,1]])
sage: M2 = CoxeterMatrix([[1,3,3,3],[3,1,3,3],[3,3,1,3],[3,3,3,1]])
sage: M3 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,3,3],[oo,3,1,3],[oo,3,3,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: GR2 = GeometricRepresentationCoxeterGroup(M2)
sage: GR3 = GeometricRepresentationCoxeterGroup(M3)
sage: GR1.fundamental_space_weights()
[(1/4, -1/4, -1/4, -1/4),
(-1/4, 1/4, -1/4, -1/4),
(-1/4, -1/4, 1/4, -1/4),
(-1/4, -1/4, -1/4, 1/4)]
sage: GR2.fundamental_space_weights()
[]
sage: GR3.fundamental_space_weights()
[(-1/3, 1/3, -1/3, -1/3), (-1/3, -1/3, 1/3, -1/3), (-1/3, -1/3, -1/3, 1/3)]
"""
fund_weights = self.fundamental_weights(exact)
bf = self.bilinear_form(exact)
fundamental_space_weights = [fw for fw in fund_weights if fw*bf*fw > 0]
self._computed_space_weights = Set(fundamental_space_weights)
return fundamental_space_weights
[docs] @cached_method
def signature(self):
r"""
Return the signature of the associated quadratic space, i.e., the vector space equipped
with the bilinear form
OUTPUT:
A 3-tuple containing the number of positive, negative and zero
eigenvalues of the bilinear form
EXAMPLES:
A finite type::
sage: from brocoli import *
sage: B2 = CoxeterMatrix([[1,4],[4,1]])
sage: GR_B2 = GeometricRepresentationCoxeterGroup(B2)
sage: GR_B2.signature()
(2, 0, 0)
An affine type has zero as an eigenvalue of multiplicity 1::
sage: Di_Af = CoxeterMatrix([[1,oo],[oo,1]])
sage: GR_DiAf = GeometricRepresentationCoxeterGroup(Di_Af)
sage: GR_DiAf.signature()
(1, 0, 1)
A Lorentzian type in rank 2::
sage: Di_Hy = CoxeterMatrix([[1,-2],[-2,1]])
sage: GR_DiHy = GeometricRepresentationCoxeterGroup(Di_Hy)
sage: GR_DiHy.signature()
(1, 1, 0)
A Lorentzian type in rank 4::
sage: M1 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,oo,oo],[oo,oo,1,oo],[oo,oo,oo,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: GR1.signature()
(3, 1, 0)
The affine type `\tilde{C_4}`::
sage: C4t = CoxeterMatrix([[1,4,2,2],[4,1,3,2],[2,3,1,4],[2,2,4,1]])
sage: GR_C4t = GeometricRepresentationCoxeterGroup(C4t)
sage: GR_C4t.signature()
(3, 0, 1)
A special Lorentzian type::
sage: L1 = CoxeterMatrix([[1,oo,2,2],[oo,1,3,2],[2,3,1,oo],[2,2,oo,1]])
sage: GR_L1 = GeometricRepresentationCoxeterGroup(L1)
sage: GR_L1.signature()
(3, 1, 0)
If we change the entry for the infinite labels to `-\sqrt{3}/2`, we
get::
sage: QF.<a> = QuadraticField(3/2)
sage: Cyl = CoxeterMatrix([[1,-a,2,2],[-a,1,3,2],[2,3,1,-a],[2,2,-a,1]])
sage: GR_Cyl = GeometricRepresentationCoxeterGroup(Cyl)
sage: GR_Cyl.signature()
(2, 1, 1)
Putting a small label gives::
sage: L2 = CoxeterMatrix([[1,-3/2,2,2],[-3/2,1,3,2],[2,3,1,-3/2],[2,2,-3/2,1]])
sage: GR_L2 = GeometricRepresentationCoxeterGroup(L2)
sage: GR_L2.signature()
(2, 2, 0)
Using an inexact base ring raises an error::
sage: L3 = CoxeterMatrix([[1,-1.5,2,2],[-1.5,1,3,2],[2,3,1,-1.5],[2,2,-1.5,1]])
sage: GR_L3 = GeometricRepresentationCoxeterGroup(L3)
sage: GR_L3.signature()
Traceback (most recent call last):
...
ValueError: the base ring is not exact
TESTS::
sage: Finite = [GeometricRepresentationCoxeterGroup(CM) for CM in CoxeterMatrix.samples(finite=True)]
sage: all(gm.signature()[1:] == (0,0) for gm in Finite)
True
sage: Affine = [GeometricRepresentationCoxeterGroup(CM) for CM in CoxeterMatrix.samples(affine=True)]
sage: all(gm.signature()[1:] == (0,1) for gm in Affine)
True
"""
if self._base_ring is RDF:
raise ValueError("the base ring is not exact")
else:
eigen_bilin_form = (self.bilinear_form()).change_ring(QQbar)
ev_signs = [sign(m.real()) for m in eigen_bilin_form.eigenvalues()]
return (ev_signs.count(1), ev_signs.count(-1), ev_signs.count(0))
[docs] def is_finite(self):
"""
Check if ``self`` is a geometric representation of a finite Coxeter group
EXAMPLES::
sage: from brocoli import *
sage: B2 = CoxeterMatrix([[1,4],[4,1]])
sage: A3 = CoxeterMatrix([[1,3,2],[3,1,3],[2,3,1]])
sage: GR_B2 = GeometricRepresentationCoxeterGroup(B2)
sage: GR_A3 = GeometricRepresentationCoxeterGroup(A3)
sage: GR_B2.is_finite()
True
sage: GR_A3.is_finite()
True
sage: Di_Af = CoxeterMatrix([[1,oo],[oo,1]])
sage: GR_DiAf = GeometricRepresentationCoxeterGroup(Di_Af)
sage: GR_DiAf.is_finite()
False
TESTS::
sage: Finite = [GeometricRepresentationCoxeterGroup(CM) for CM in CoxeterMatrix.samples(finite=True)]
sage: all([f.is_finite() for f in Finite])
True
sage: Affine = [GeometricRepresentationCoxeterGroup(CM) for CM in CoxeterMatrix.samples(affine=True)]
sage: any([a.is_finite() for a in Affine])
False
"""
sig = self.signature()
if sig[1] == 0 and sig[2] == 0: # Finite type
return True
else:
return False
[docs] def is_affine(self):
"""
Check if ``self`` is a geometric representation of an affine Coxeter group
EXAMPLES::
sage: from brocoli import *
sage: C4t = CoxeterMatrix([[1,4,2,2],[4,1,3,2],[2,3,1,4],[2,2,4,1]])
sage: GR_C4t = GeometricRepresentationCoxeterGroup(C4t)
sage: GR_C4t.is_affine()
True
sage: Di_Af = CoxeterMatrix([[1,oo],[oo,1]])
sage: GR_DiAf = GeometricRepresentationCoxeterGroup(Di_Af)
sage: GR_DiAf.is_affine()
True
sage: B2 = CoxeterMatrix([[1,4],[4,1]])
sage: GR_B2 = GeometricRepresentationCoxeterGroup(B2)
sage: GR_B2.is_affine()
False
TESTS::
sage: from brocoli import *
sage: Affine = [GeometricRepresentationCoxeterGroup(CM) for CM in CoxeterMatrix.samples(affine=True)]
sage: all([a.is_affine() for a in Affine])
True
sage: Finite = [GeometricRepresentationCoxeterGroup(CM) for CM in CoxeterMatrix.samples(finite=True)]
sage: any([f.is_affine() for f in Finite])
False
"""
sig = self.signature()
if sig[1] == 0 and sig[2] == 1: # Affine type
return True
else:
return False
[docs] def is_lorentzian(self):
r"""
Check if ``self`` is a geometric representation of a Lorentzian Coxeter group
EXAMPLES::
sage: from brocoli import *
sage: L1 = CoxeterMatrix([[1,oo,2,2],[oo,1,3,2],[2,3,1,oo],[2,2,oo,1]])
sage: L2 = CoxeterMatrix([[1,-3/2,2,2],[-3/2,1,3,2],[2,3,1,-3/2],[2,2,-3/2,1]])
sage: GR_L1 = GeometricRepresentationCoxeterGroup(L1)
sage: GR_L2 = GeometricRepresentationCoxeterGroup(L2)
sage: GR_L1.is_lorentzian()
True
sage: GR_L2.is_lorentzian()
False
If the bilinear form has only one negative eigenvalue, it is not
necessarily lorentzian::
sage: QF.<a> = QuadraticField(3/2)
sage: Cyl = CoxeterMatrix([[1,-a,2,2],[-a,1,3,2],[2,3,1,-a],[2,2,-a,1]])
sage: GR_Cyl = GeometricRepresentationCoxeterGroup(Cyl)
sage: GR_Cyl.is_lorentzian()
False
sage: GR_Cyl.signature()
(2, 1, 1)
.. SEEALSO::
- :meth:`is_degenerate_lorentzian`
TESTS::
sage: Finite = [GeometricRepresentationCoxeterGroup(CM) for CM in CoxeterMatrix.samples(finite=True)]
sage: all([not a.is_lorentzian() for a in Finite])
True
sage: Affine = [GeometricRepresentationCoxeterGroup(CM) for CM in CoxeterMatrix.samples(affine=True)]
sage: all([not a.is_lorentzian() for a in Affine])
True
"""
sig = self.signature()
if sig[1] == 1 and sig[2] == 0: # Lorentzian type
return True
else:
return False
[docs] def is_degenerate_lorentzian(self):
r"""
Check if ``self`` is a geometric representation of a Lorentzian Coxeter group with a
singular bilinear form
EXAMPLES::
sage: from brocoli import *
sage: QF.<a> = QuadraticField(3/2)
sage: Cyl = CoxeterMatrix([[1,-a,2,2],[-a,1,3,2],[2,3,1,-a],[2,2,-a,1]])
sage: GR_Cyl = GeometricRepresentationCoxeterGroup(Cyl)
sage: GR_Cyl.is_degenerate_lorentzian()
True
"""
sig = self.signature()
if sig[1] == 1 and sig[2] > 0: # Singular Lorentzian type
return True
else:
return False
[docs] @cached_method
def elements(self, length):
r"""
Return the elements of the represented Coxeter group with length given by ``length``.
INPUT:
- ``length`` -- a non-negative integer
OUTPUT:
A set of matrices
EXAMPLES:
Finite examples::
sage: from brocoli import *
sage: Di5 = CoxeterMatrix([[1,5],[5,1]])
sage: GR_Di5 = GeometricRepresentationCoxeterGroup(Di5)
sage: [len(GR_Di5.elements(i)) for i in range(7)]
[1, 2, 2, 2, 2, 1, 0]
sage: A3 = CoxeterMatrix([[1,3,2],[3,1,3],[2,3,1]])
sage: GR_A3 = GeometricRepresentationCoxeterGroup(A3)
sage: [len(GR_A3.elements(i)) for i in range(8)]
[1, 3, 5, 6, 5, 3, 1, 0]
Infinite examples::
sage: M1 = CoxeterMatrix([[1,4,4],[4,1,4],[4,4,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: [len(GR1.elements(i)) for i in range(6)]
[1, 3, 6, 12, 21, 36]
sage: M2 = CoxeterMatrix([[1,4,-3/2],[4,1,4],[-3/2,4,1]])
sage: GR2 = GeometricRepresentationCoxeterGroup(M2)
sage: [len(GR2.elements(i)) for i in range(6)]
[1, 3, 6, 12, 22, 40]
TESTS::
sage: GR1.elements(-1)
Traceback (most recent call last):
...
ValueError: length has to be a positive integer
"""
if length < 0:
raise ValueError("length has to be a positive integer")
elif length == 0:
return Set([self.identity_element()])
elif length == 1:
return Set(self.simple_reflections())
else:
element_set = Set([])
for element in self.elements(length - 1):
for simple_refl in self.simple_reflections():
New_element = simple_refl*element
New_element.set_immutable()
if New_element not in self.elements(length - 2):
new_word = self._matrix_to_word[element] + self._matrix_to_word[simple_refl]
self._matrix_to_word[New_element] = new_word
element_set = element_set.union(Set([New_element]))
return element_set
@cached_method
def _real_element(self, element):
r"""
Return a floating point approximation of the ``element``
INPUT:
- ``element`` -- matrix; a matrix representing an element of the representation
OUTPUT:
a matrix over ``RDF``
EXAMPLES::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,4,4],[4,1,4],[4,4,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: GR._real_element(GR.identity_element())
[1.0 0.0 0.0]
[0.0 1.0 0.0]
[0.0 0.0 1.0]
sage: GR._real_element(GR.elements(2)[3])
[ -1.0 1.414213562373095 1.414213562373095]
[ 0.0 1.0 0.0]
[-1.414213562373095 3.414213562373095 1.0]
"""
rdf_elmt = element.change_ring(RDF)
rdf_elmt.set_immutable()
return rdf_elmt
@cached_method
def _algebraic_element(self, element):
r"""
Return ``element`` in the base ring ``QQbar``.
INPUT:
- ``element`` -- matrix; a matrix representing an element of the representation
OUTPUT:
a matrix over ``QQbar``
EXAMPLES::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,4,4],[4,1,4],[4,4,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: GR._algebraic_element(GR.identity_element()).base_ring()
Algebraic Field
sage: GR._algebraic_element(GR.elements(2)[4])
[ 1 0 0]
[ 3.414213562373095? + 0.?e-18*I 1 -1.414213562373095? + 0.?e-18*I]
[ 1.414213562373095? + 0.?e-18*I 1.414213562373095? + 0.?e-18*I -1]
.. NOTE::
This is used to compute minimal polynomial and eigenvalues. The
minimal polynomial can not be computed over the Cyclotomic Field
yet.
"""
if not self._base_ring.is_exact():
exact_elmt = element.change_ring(QQbar)
exact_elmt.set_immutable()
return exact_elmt
elif self._base_ring is UniversalCyclotomicField(): # TODO: Get rid of this hack
exact_elmt = element.change_ring(QQbar)
exact_elmt.set_immutable()
return exact_elmt
else:
return element
[docs] def matrix_to_word(self, matrix):
r"""
Return a reduced word associated to ``matrix``
EXAMPLES::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,oo,oo,oo],[oo,1,oo,oo],[oo,oo,1,oo],[oo,oo,oo,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
The matrices have to be created beforehand inside the representation::
sage: my_ident = identity_matrix(GR.base_ring(),4)
sage: GR.matrix_to_word(my_ident)
Traceback (most recent call last):
...
ValueError: the matrix was not computed before.
sage: e = GR.identity_element()
sage: GR.matrix_to_word(e)
word:
sage: elmts = GR.elements(2)
sage: for elmt in elmts:
....: print(GR.matrix_to_word(elmt))
....:
24
41
14
43
13
23
34
42
21
31
12
32
"""
if matrix in self._matrix_to_word.keys():
return self._matrix_to_word[matrix]
else:
raise ValueError("the matrix was not computed before.")
@cached_method
def _minimal_polynomial(self, element):
r"""
Return the minimal polynomial of ``element``
INPUT:
- ``element`` -- matrix; a matrix representing an element of the representation
OUTPUT:
a polynomial
EXAMPLES::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,3,oo,oo],[3,1,oo,oo],[oo,oo,1,oo],[oo,oo,oo,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: e = GR.identity_element()
sage: elliptic = GR.simple_reflections()[0]*GR.simple_reflections()[1]
sage: parabolic = GR.simple_reflections()[2]*GR.simple_reflections()[3]
sage: hyperbolic1 = GR.simple_reflections()[1]*GR.simple_reflections()[2]*GR.simple_reflections()[3]
sage: hyperbolic2 = GR.elements(4)[0]
sage: GR._minimal_polynomial(elliptic)
x^3 - 1
sage: GR._minimal_polynomial(parabolic)
x^3 - 3*x^2 + 3*x - 1
sage: GR._minimal_polynomial(hyperbolic1)
x^4 - 18*x^3 + 18*x - 1
sage: GR._minimal_polynomial(hyperbolic2)
x^4 - 49*x^3 - 96*x^2 - 49*x + 1
"""
return self._algebraic_element(element).minimal_polynomial()
@cached_method
def _is_diagonalizable(self, element):
r"""
Check if ``element`` is diagonalizable. This is equivalent to the fact that its
minimal polynomial is squarefree
INPUT:
- ``element`` -- a matrix coming from the representation
EXAMPLES::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,3,oo,oo],[3,1,oo,oo],[oo,oo,1,oo],[oo,oo,oo,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: e = GR.identity_element()
sage: elliptic = GR.simple_reflections()[0]*GR.simple_reflections()[1]
sage: parabolic = GR.simple_reflections()[2]*GR.simple_reflections()[3]
sage: hyperbolic1 = GR.simple_reflections()[1]*GR.simple_reflections()[2]*GR.simple_reflections()[3]
sage: hyperbolic2 = GR.elements(4)[0]
sage: GR._is_diagonalizable(e)
True
sage: GR._is_diagonalizable(elliptic)
True
sage: GR._is_diagonalizable(parabolic)
False
sage: GR._is_diagonalizable(hyperbolic1)
True
sage: GR._is_diagonalizable(hyperbolic2)
True
"""
return self._minimal_polynomial(element).is_squarefree()
@cached_method
def _algebraic_eigenvalues(self, element):
r"""
Return the eigenvalues of ``element``
INPUT:
- ``element`` -- matrix; a matrix representing an element of the representation
OUTPUT:
A list of pairs `(r,m)` giving the roots and their multiplicities
EXAMPLES::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,3,oo,oo],[3,1,oo,oo],[oo,oo,1,oo],[oo,oo,oo,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: e = GR.identity_element()
sage: elliptic = GR.simple_reflections()[0]*GR.simple_reflections()[1]
sage: parabolic = GR.simple_reflections()[2]*GR.simple_reflections()[3]
sage: hyperbolic1 = GR.simple_reflections()[1]*GR.simple_reflections()[2]*GR.simple_reflections()[3]
sage: hyperbolic2 = GR.elements(4)[0]
sage: GR._algebraic_eigenvalues(e)
[(1, 1)]
sage: GR._algebraic_eigenvalues(elliptic)
[(1, 1),
(-0.500000000000000? - 0.866025403784439?*I, 1),
(-0.500000000000000? + 0.866025403784439?*I, 1)]
sage: GR._algebraic_eigenvalues(parabolic)
[(1, 3)]
sage: GR._algebraic_eigenvalues(hyperbolic1)
[(-1, 1), (0.05572809000084122?, 1), (1, 1), (17.94427190999916?, 1)]
sage: GR._algebraic_eigenvalues(hyperbolic2)
[(0.01964452215592945?, 1),
(50.90477600129170?, 1),
(-0.9622102617238120? - 0.2723075691812354?*I, 1),
(-0.9622102617238120? + 0.2723075691812354?*I, 1)]
sage: M1 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,3,3],[oo,3,1,3],[oo,3,3,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: a,b,c,d = GR1.simple_reflections()
sage: my_elmt = c*d*a*b
sage: GR._algebraic_eigenvalues(my_elmt)
[(0.03351094477932432?, 1),
(29.84099692160822?, 1),
(-0.9372539331937718? - 0.3486474791430519?*I, 1),
(-0.9372539331937718? + 0.3486474791430519?*I, 1)]
"""
return self._minimal_polynomial(element).roots(ring=QQbar)
@cached_method
def _algebraic_eigenvectors_element(self, element):
"""
Return the eigenvectors of ``element`` in an exact ring
INPUT:
- ``element`` -- matrix; a matrix representing an element of the representation
OUTPUT:
eigenvectors of ``element``
EXAMPLES::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,4,oo,oo],[4,1,oo,oo],[oo,oo,1,oo],[oo,oo,oo,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: e = GR.identity_element()
sage: elliptic = GR.simple_reflections()[0]*GR.simple_reflections()[1]
sage: parabolic = GR.simple_reflections()[2]*GR.simple_reflections()[3]
sage: hyperbolic = GR.simple_reflections()[1]*GR.simple_reflections()[2]*GR.simple_reflections()[3]
sage: GR._algebraic_eigenvectors_element(e)
[(1, [
(1, 0, 0, 0),
(0, 1, 0, 0),
(0, 0, 1, 0),
(0, 0, 0, 1)
], 4)]
sage: GR._algebraic_eigenvectors_element(elliptic)
[(1*I, [
(1, 0.7071067811865475? - 0.7071067811865475?*I, 0, 0)
], 1), (-1*I, [
(1, 0.7071067811865475? + 0.7071067811865475?*I, 0, 0)
], 1), (1, [
(1, 1, 0, 0.2928932188134525? + 0.?e-18*I),
(0, 0, 1, -1.000000000000000? + 0.?e-18*I)
], 2)]
sage: GR._algebraic_eigenvectors_element(parabolic)
[(1, [
(1, -1, 0, 0),
(0, 0, 1.000000000000000?, 1.000000000000000?)
], 4)]
sage: GR._algebraic_eigenvectors_element(hyperbolic)
[(17.94427190999916?, [
(0, 1.000000000000000?, 0.381966011250106?, 0.1458980337503155?)
], 1), (1, [
(1.000000000000000? + 0.?e-17*I, -1.000000000000000? + 0.?e-17*I, -0.8535533905932738? + 0.?e-17*I, -0.8535533905932738? + 0.?e-17*I)
], 1), (0.05572809000084122?, [
(0, 1.000000000000000?, 2.618033988749895?, 6.854101966249684?)
], 1), (-1, [
(0, 1, -1, 1)
], 1)]
"""
M = self._algebraic_element(element)
return M.eigenvectors_right()
[docs] @cached_method
def roots(self, depth):
r"""
Return the roots of the representation of depth ``depth``
INPUT:
- ``depth`` -- a non-negative integer
OUTPUT:
A set of roots of depth ``depth``
EXAMPLES::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,4,oo,oo],[4,1,oo,oo],[oo,oo,1,-2],[oo,oo,-2,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: GR.roots(0)
{(0, 0, 0, 1), (1, 0, 0, 0), (0, 1, 0, 0), (0, 0, 1, 0)}
sage: len(GR.roots(1))
12
sage: len(GR.roots(2))
34
sage: A3 = CoxeterMatrix([[1,3,2],[3,1,3],[2,3,1]])
sage: GR_A3 = GeometricRepresentationCoxeterGroup(A3)
sage: [len(GR_A3.roots(i)) for i in range(4)]
[3, 2, 1, 0]
sage: QF.<a> = QuadraticField(3/2)
sage: M2 = CoxeterMatrix([[1,-a,5],[-a,1,-3/2],[5,-3/2,1]])
sage: GR2 = GeometricRepresentationCoxeterGroup(M2)
sage: [len(GR2.roots(i)) for i in range(5)]
[3, 4, 7, 11, 18]
"""
if depth < 0:
raise ValueError("depth has to be a positive integer.")
elif depth == 0:
return Set(self.simple_roots())
else:
self.roots(depth-1)
set_roots = Set([])
for element in self.elements(depth):
for col in range(self._rank):
nr = element.column(col)
if nr < 0:
nr = -nr
nr.set_immutable()
if nr not in self._computed_roots and nr not in set_roots:
set_roots = set_roots.union(Set([nr]))
self._computed_roots = self._computed_roots.union(set_roots)
return set_roots
@cached_method
def _rdf_vector(self, v):
r"""
Return a floating point approximation of the vector``v``
INPUT:
- ``v`` -- vector;
OUTPUT:
A vector over ``RDF``
EXAMPLES::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,4,oo,oo],[4,1,oo,oo],[oo,oo,1,-2],[oo,oo,-2,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: roots = GR.roots(2)
sage: for r in roots:
....: print(GR._rdf_vector(r)) # abs tol 1e-14
....:
(5.414213562373096, 1.0, 0.0, 2.0)
(0.0, 3.0, 2.0, 0.0)
(1.0, 1.414213562373095, 4.82842712474619, 0.0)
(0.0, 10.0, 4.0, 1.0)
(2.0, 4.82842712474619, 1.0, 0.0)
(4.82842712474619, 2.0, 0.0, 1.0)
(2.0, 0.0, 3.0, 0.0)
(1.0, 0.0, 10.0, 2.0)
(0.0, 2.0, 8.0, 1.0)
(1.414213562373095, 1.0, 4.82842712474619, 0.0)
(0.0, 2.0, 3.0, 0.0)
(2.0, 0.0, 1.0, 8.0)
(1.0, 0.0, 2.0, 10.0)
(2.0, 4.82842712474619, 0.0, 1.0)
(0.0, 1.0, 10.0, 2.0)
(1.0, 5.414213562373096, 0.0, 2.0)
(5.414213562373096, 1.0, 2.0, 0.0)
(1.414213562373095, 1.0, 0.0, 4.82842712474619)
(0.0, 0.0, 4.0, 15.0)
(0.0, 2.0, 0.0, 3.0)
(0.0, 0.0, 15.0, 4.0)
(0.0, 3.0, 0.0, 2.0)
(3.0, 0.0, 2.0, 0.0)
(1.0, 5.414213562373096, 2.0, 0.0)
(0.0, 1.0, 2.0, 10.0)
(10.0, 0.0, 1.0, 4.0)
(2.0, 0.0, 0.0, 3.0)
(0.0, 2.0, 1.0, 8.0)
(0.0, 10.0, 1.0, 4.0)
(4.82842712474619, 2.0, 1.0, 0.0)
(1.0, 1.414213562373095, 0.0, 4.82842712474619)
(10.0, 0.0, 4.0, 1.0)
(2.0, 0.0, 8.0, 1.0)
(3.0, 0.0, 0.0, 2.0)
sage: [GR._rdf_vector(w) for w in GR.weights(2)] # abs tol 1e-14
[(-0.23294313392598154, -0.23294313392598154, -8.200780398814702, -1.8674470654813702),
(-0.1764216518504617, -1.590635214223557, -4.232943133925982, -0.23294313392598154),
(-0.23294313392598154, -0.23294313392598154, -1.8674470654813702, -8.200780398814702),
(-3.004848776596652, -1.590635214223557, -0.23294313392598154, -0.23294313392598154),
(-0.1764216518504617, -1.590635214223557, -0.23294313392598154, -4.232943133925982),
(-1.590635214223557, -3.004848776596652, -0.23294313392598154, -0.23294313392598154),
(-0.23294313392598154, -4.232943133925982, -1.8674470654813702, -0.2007803988147036),
(-1.590635214223557, -0.1764216518504617, -0.23294313392598154, -4.232943133925982),
(-4.232943133925982, -0.23294313392598154, -1.8674470654813702, -0.2007803988147036),
(-4.232943133925982, -0.23294313392598154, -0.2007803988147036, -1.8674470654813702),
(-0.23294313392598154, -4.232943133925982, -0.2007803988147036, -1.8674470654813702),
(-1.590635214223557, -0.1764216518504617, -4.232943133925982, -0.23294313392598154)]
sage: l_roots = GR.limit_roots(2) # long time
sage: for lr in l_roots: # long time
....: print(GR._rdf_vector(lr)) # long time abs tol 1e-14
....:
(0.0, 1.0, 1.0, -0.0)
(0.0, 1.0, -0.0, 1.0)
(1.0, 0.0, 1.0, -0.0)
(0.0, 0.0, 1.0, 3.7320508075688776)
(1.0, 0.0, -0.0, 1.0)
(0.0, 0.0, 1.0, 0.2679491924311227)
"""
rdf_vector = v.change_ring(RDF)
rdf_vector.set_immutable()
return rdf_vector
[docs] @cached_method
def weights(self, depth, exact=True):
r"""
Return the weights of the representation of depth ``depth``
INPUT:
- ``depth`` -- a non-negative integer
- ``exact`` -- a boolean (default = ``True``); whether to give the bilinear in an exact base
ring or not.
OUTPUT:
A set of weights of depth ``depth``
EXAMPLES::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,4,oo,oo],[4,1,oo,oo],[oo,oo,1,-2],[oo,oo,-2,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: [len(GR.weights(i)) for i in range(3)]
[4, 4, 12]
Putting ``exact=False`` will give more weights as computations are not
exact. This may be useful if we just want to draw weights::
sage: [len(GR.weights(i,False)) for i in range(3)] # not tested
[4, 16, 48]
A finite example::
sage: A3 = CoxeterMatrix([[1,3,2],[3,1,3],[2,3,1]])
sage: GR_A3 = GeometricRepresentationCoxeterGroup(A3)
sage: [len(GR_A3.weights(i)) for i in range(6)]
[3, 3, 4, 3, 1, 0]
"""
if depth < 0:
raise ValueError("depth has to be a positive integer.")
elif depth == 0:
return Set(self.fundamental_weights(exact))
else:
self.weights(depth - 1)
if exact:
fw_matrix = self._fund_weights_matrix_exact
else:
fw_matrix = self._fund_weights_matrix_rdf
set_weights = Set([])
for element in self.elements(depth):
if exact:
weight_basis = element*fw_matrix
else:
weight_basis = self._real_element(element)*fw_matrix
for col in range(self._rank):
new_weight = weight_basis.column(col)
new_weight.set_immutable()
if exact:
if new_weight not in self._computed_weights_exact and new_weight not in set_weights:
set_weights = set_weights.union(Set([new_weight]))
else:
if new_weight not in self._computed_weights_rdf and new_weight not in set_weights:
set_weights = set_weights.union(Set([new_weight]))
if exact:
self._computed_weights_exact = self._computed_weights_exact.union(set_weights)
else:
self._computed_weights_rdf = self._computed_weights_rdf.union(set_weights)
return set_weights
[docs] @cached_method
def space_weights(self, depth, exact=True):
r"""
Return the weights of the representation with positive bilinear form
evaluation
INPUT:
- ``depth`` -- a non-negative integer
- ``exact`` -- a boolean (default = ``True``); whether to give the bilinear in an exact base
ring or not.
OUTPUT:
A set of weights
EXAMPLES::
sage: from brocoli import *
sage: M1 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,oo,oo],[oo,oo,1,oo],[oo,oo,oo,1]])
sage: M2 = CoxeterMatrix([[1,3,3,3],[3,1,3,3],[3,3,1,3],[3,3,3,1]])
sage: M3 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,3,3],[oo,3,1,3],[oo,3,3,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: GR2 = GeometricRepresentationCoxeterGroup(M2)
sage: GR3 = GeometricRepresentationCoxeterGroup(M3)
sage: GR1.space_weights(0)
{(-1/4, 1/4, -1/4, -1/4), (-1/4, -1/4, 1/4, -1/4), (1/4, -1/4, -1/4, -1/4), (-1/4, -1/4, -1/4, 1/4)}
sage: GR1.space_weights(1)
{(-1/4, -7/4, -1/4, -1/4), (-1/4, -1/4, -7/4, -1/4), (-7/4, -1/4, -1/4, -1/4), (-1/4, -1/4, -1/4, -7/4)}
sage: [len(GR1.weights(i)) for i in range(5)]
[4, 4, 12, 36, 108]
sage: GR2.space_weights(0)
{}
sage: GR3.space_weights(0)
{(-1/3, 1/3, -1/3, -1/3), (-1/3, -1/3, -1/3, 1/3), (-1/3, -1/3, 1/3, -1/3)}
sage: [len(GR3.weights(i)) for i in range(5)]
[4, 4, 12, 30, 84]
sage: [len(GR3.weights(i,False)) for i in range(5)] # not tested
[4, 8, 30, 91, 273]
"""
if depth < 0:
raise ValueError("depth has to be a positive integer.")
elif depth == 0:
return Set(self.fundamental_space_weights(exact))
else:
set_space_weights = Set([fw for fw in self.weights(depth, exact) if fw*self.bilinear_form(exact)*fw > 0])
if exact:
self._computed_space_weights_exact = self._computed_space_weights_exact.union(set_space_weights)
else:
self._computed_space_weights_rdf = self._computed_space_weights_rdf.union(set_space_weights)
return set_space_weights
@cached_method
def _split_order_elements(self, length):
r"""
For Lorentzian representation, partitions the elements of length
``length`` into the three type of transformations of Lorentz space:
- Elliptic: Finite order elements
- Parabolic: Infinite but not diagonalizable
- Hyperbolic: Infinite and diagonalizable
INPUT:
- ``length`` -- a non-negative integer
OUTPUT:
A triple of sets containing the elliptic, parabolic and hyperbolic
elements respectively.
EXAMPLES::
sage: from brocoli import *
sage: M1 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,oo,oo],[oo,oo,1,oo],[oo,oo,oo,1]])
sage: M2 = CoxeterMatrix([[1,3,3,3],[3,1,3,3],[3,3,1,3],[3,3,3,1]])
sage: M3 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,3,3],[oo,3,1,3],[oo,3,3,1]])
sage: M4 = CoxeterMatrix([[1,-2,-2,-2],[-2,1,-2,-2],[-2,-2,1,-2],[-2,-2,-2,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: GR2 = GeometricRepresentationCoxeterGroup(M2)
sage: GR3 = GeometricRepresentationCoxeterGroup(M3)
sage: GR4 = GeometricRepresentationCoxeterGroup(M4)
sage: for gr in [GR1,GR2,GR3,GR4]:
....: for ell in range(4):
....: split = gr._split_order_elements(ell)
....: print(ell,[len(s) for s in split])
....:
0 [1, 0, 0]
1 [4, 0, 0]
2 [0, 12, 0]
3 [12, 0, 24]
0 [1, 0, 0]
1 [4, 0, 0]
2 [12, 0, 0]
3 [6, 24, 0]
0 [1, 0, 0]
1 [4, 0, 0]
2 [6, 6, 0]
3 [9, 6, 18]
0 [1, 0, 0]
1 [4, 0, 0]
2 [0, 0, 12]
3 [12, 0, 24]
sage: M5 = CoxeterMatrix([[1,4,oo,oo],[4,1,oo,oo],[oo,oo,1,-2],[oo,oo,-2,1]])
sage: GR5 = GeometricRepresentationCoxeterGroup(M5)
sage: [len(i) for i in GR5._split_order_elements(3)]
[12, 0, 24]
sage: QF.<a> = QuadraticField(3/2)
sage: Cyl = CoxeterMatrix([[1,-a,2,2],[-a,1,3,2],[2,3,1,-a],[2,2,-a,1]])
sage: GR_Cyl = GeometricRepresentationCoxeterGroup(Cyl)
sage: GR_Cyl._split_order_elements(0)
Traceback (most recent call last):
...
NotImplementedError: The representation space has to be Euclidean, Affine or Lorentzian to split the elements
sage: GR1._split_order_elements(-1)
Traceback (most recent call last):
...
ValueError: length has to be a positive integer >=0
"""
if not self.is_finite() and not self.is_affine() and not self.is_lorentzian():
raise NotImplementedError("The representation space has to be Euclidean, Affine or Lorentzian to split the elements")
elif length < 0:
raise ValueError("length has to be a positive integer >=0")
elif length == 0:
return self.elements(0), Set([]), Set([])
else:
elliptic_elements = Set([])
parabolic_elements = Set([])
hyperbolic_elements = Set([])
for element in self.elements(length):
if self._is_diagonalizable(element):
Eigenvalues = self._algebraic_eigenvalues(element)
if max([ev[0].norm() for ev in Eigenvalues]) > 1.01: # Element is hyperbolic
hyperbolic_elements = hyperbolic_elements.union(Set([element]))
else: # Element is elliptic
elliptic_elements = elliptic_elements.union(Set([element]))
else: # Element is parabolic
parabolic_elements = parabolic_elements.union(Set([element]))
return elliptic_elements, parabolic_elements, hyperbolic_elements
[docs] def elliptic_elements(self, length):
r"""
For Lorentzian Coxeter groups, return the elliptic elements of the
representation of length ``length``
INPUT:
- ``length`` -- a non-negative integer
OUTPUT:
A set of elements
EXAMPLES::
sage: from brocoli import *
sage: M1 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,oo,oo],[oo,oo,1,oo],[oo,oo,oo,1]])
sage: M2 = CoxeterMatrix([[1,3,3,3],[3,1,3,3],[3,3,1,3],[3,3,3,1]])
sage: M3 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,3,3],[oo,3,1,3],[oo,3,3,1]])
sage: M4 = CoxeterMatrix([[1,-2,-2,-2],[-2,1,-2,-2],[-2,-2,1,-2],[-2,-2,-2,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: GR2 = GeometricRepresentationCoxeterGroup(M2)
sage: GR3 = GeometricRepresentationCoxeterGroup(M3)
sage: GR4 = GeometricRepresentationCoxeterGroup(M4)
sage: GR1.elliptic_elements(0)
{[1 0 0 0]
[0 1 0 0]
[0 0 1 0]
[0 0 0 1]}
sage: len(GR1.elliptic_elements(3))
12
sage: len(GR1.elliptic_elements(4))
0
sage: len(GR2.elliptic_elements(3))
6
sage: len(GR2.elliptic_elements(4))
24
sage: len(GR3.elliptic_elements(3))
9
sage: len(GR3.elliptic_elements(4))
12
sage: GR4.elliptic_elements(4)[-1]
Traceback (most recent call last):
...
IndexError: list index out of range
TESTS::
sage: for gr in [GR1,GR2,GR3,GR4]:
....: for ell in [0,4]:
....: elliptic = gr.elliptic_elements(ell)
....: print(ell, len(elliptic))
....: if len(elliptic)>0:
....: print(elliptic[-1])
....:
0 1
[1 0 0 0]
[0 1 0 0]
[0 0 1 0]
[0 0 0 1]
4 0
0 1
[1 0 0 0]
[0 1 0 0]
[0 0 1 0]
[0 0 0 1]
4 24
[ 0 -1 2 2]
[ 0 -2 6 3]
[ 0 0 1 0]
[ 1 -2 4 2]
0 1
[1 0 0 0]
[0 1 0 0]
[0 0 1 0]
[0 0 0 1]
4 12
[ 1 0 0 0]
[12 -2 0 3]
[ 4 -1 0 2]
[ 8 -2 1 2]
0 1
[1 0 0 0]
[0 1 0 0]
[0 0 1 0]
[0 0 0 1]
4 0
"""
return self._split_order_elements(length)[0]
[docs] def parabolic_elements(self, length):
r"""
For Lorentzian Coxeter groups, return the parabolic elements of the
representation of length ``length``
INPUT:
- ``length`` -- a non-negative integer
OUTPUT:
A set of elements
EXAMPLES::
sage: from brocoli import *
sage: M1 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,oo,oo],[oo,oo,1,oo],[oo,oo,oo,1]])
sage: M2 = CoxeterMatrix([[1,3,3,3],[3,1,3,3],[3,3,1,3],[3,3,3,1]])
sage: M3 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,3,3],[oo,3,1,3],[oo,3,3,1]])
sage: M4 = CoxeterMatrix([[1,-2,-2,-2],[-2,1,-2,-2],[-2,-2,1,-2],[-2,-2,-2,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: GR2 = GeometricRepresentationCoxeterGroup(M2)
sage: GR3 = GeometricRepresentationCoxeterGroup(M3)
sage: GR4 = GeometricRepresentationCoxeterGroup(M4)
sage: GR1.parabolic_elements(0)
{}
sage: len(GR1.parabolic_elements(3))
0
sage: len(GR1.parabolic_elements(4))
36
sage: len(GR2.parabolic_elements(3))
24
sage: len(GR2.parabolic_elements(4))
24
sage: len(GR3.parabolic_elements(3))
6
sage: len(GR3.parabolic_elements(4))
24
sage: GR4.parabolic_elements(4)[-1]
Traceback (most recent call last):
...
IndexError: list index out of range
"""
return self._split_order_elements(length)[1]
[docs] def hyperbolic_elements(self, length):
r"""
For Lorentzian Coxeter groups, return the hyperbolic elements of the
representation of length ``length``
INPUT:
- ``length`` -- a non-negative integer
OUTPUT:
A set of elements
EXAMPLES::
sage: from brocoli import *
sage: M1 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,oo,oo],[oo,oo,1,oo],[oo,oo,oo,1]])
sage: M2 = CoxeterMatrix([[1,3,3,3],[3,1,3,3],[3,3,1,3],[3,3,3,1]])
sage: M3 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,3,3],[oo,3,1,3],[oo,3,3,1]])
sage: M4 = CoxeterMatrix([[1,-2,-2,-2],[-2,1,-2,-2],[-2,-2,1,-2],[-2,-2,-2,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: GR2 = GeometricRepresentationCoxeterGroup(M2)
sage: GR3 = GeometricRepresentationCoxeterGroup(M3)
sage: GR4 = GeometricRepresentationCoxeterGroup(M4)
sage: GR1.hyperbolic_elements(2)[0]
Traceback (most recent call last):
...
IndexError: list index out of range
sage: len(GR1.hyperbolic_elements(3))
24
sage: len(GR1.hyperbolic_elements(4))
72
sage: len(GR2.hyperbolic_elements(3))
0
sage: len(GR2.hyperbolic_elements(4))
24
sage: len(GR3.hyperbolic_elements(3))
18
sage: len(GR3.hyperbolic_elements(4))
54
sage: GR4.hyperbolic_elements(0)[-1]
Traceback (most recent call last):
...
IndexError: list index out of range
"""
return self._split_order_elements(length)[2]
[docs] @cached_method
def parabolic_limit_roots(self, length):
r"""
For Lorentzian Coxeter groups, return the limit roots coming from the
eigenvectors of parabolic elements of length ``length``
INPUT:
- ``length`` -- a non-negative integer
OUTPUT:
A set of vectors
EXAMPLES::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,oo,oo,oo],[oo,1,3,3],[oo,3,1,3],[oo,3,3,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: GR.parabolic_limit_roots(0)
{}
sage: GR.parabolic_limit_roots(1)
{}
sage: GR.parabolic_limit_roots(2)
{(1.0, 1.0, 0.0, -0.0), (1.0, 0.0, -0.0, 1.0), (1.0, 0.0, 1.0, -0.0)}
sage: GR.parabolic_limit_roots(3)
{(0.0, 1.0, 1.0, 1.0)}
sage: len(GR.parabolic_limit_roots(4))
10
TODO:
Make this method more robust and not return vectors of symbolic
ring
"""
para_lm = Set([])
for element in self.parabolic_elements(length):
ev_matrix = matrix(flatten([tup[1] for tup in self._algebraic_eigenvectors_element(element)])).transpose().change_ring(RDF)
if ev_matrix.ncols() == 1:
new_lr = vector(ev_matrix.column(0))
else:
extra_var = vector([var('x%i' % value) for value in range(ev_matrix.ncols())],)
ev_variables = ev_matrix*extra_var
quadric_expr = ev_variables*self.bilinear_form(False)*ev_variables
dict_soln_ev = solve([quadric_expr == 0], list(extra_var), solution_dict=True)
if len(dict_soln_ev) > 0:
new_lr = ev_variables.subs(dict_soln_ev[0])
for index in range(self._rank):
for vari in new_lr[index].variables():
new_lr[index] = new_lr[index].subs({vari: 1})
else:
soln_dict = {}
for vari in extra_var:
soln_vari = solve([quadric_expr == 0], vari)
if len(soln_vari) > 0:
soln_dict[vari] = soln_vari[0].rhs()
else:
soln_dict[vari] = 1
new_lr = ev_variables.subs(soln_dict)
if new_lr < 0:
new_lr = -new_lr
approx_problem = False
index = 0
while not approx_problem and index < self._rank:
if new_lr[index] < 0:
approx_problem = True
print "There was an approximation problem with the element", self._matrix_to_word[element]
index += 1
if not approx_problem:
new_lr.set_immutable()
if new_lr not in para_lm:
para_lm = para_lm.union(Set([new_lr]))
self._computed_limit_roots = self._computed_limit_roots.union(para_lm)
return para_lm
[docs] @cached_method
def hyperbolic_limit_roots(self, length):
r"""
For Lorentzian Coxeter groups, return the limit roots coming from the
eigenvectors of hyperbolic elements of length ``length``
INPUT:
- ``length`` -- a non-negative integer
OUTPUT:
A set of vectors
EXAMPLES::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,3,3,3],[3,1,3,3],[3,3,1,3],[3,3,3,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: [GR.hyperbolic_limit_roots(i) for i in range(4)]
[{}, {}, {}, {}]
sage: len(GR.hyperbolic_limit_roots(4)) # long time
24
sage: M1 = CoxeterMatrix([[1,-2,-2,-2],[-2,1,-2,-2],[-2,-2,1,-2],[-2,-2,-2,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: [len(GR1.hyperbolic_limit_roots(i)) for i in range(4)]
[0, 0, 12, 24]
sage: M2 = CoxeterMatrix([[1,4,4],[4,1,4],[4,4,1]])
sage: GR2 = GeometricRepresentationCoxeterGroup(M2)
sage: for lm in GR2.limit_roots(3):
....: print(lm)
....:
(1, 3.546455444684995?, 1.883203505913526?)
(1, 1.883203505913526?, 3.546455444684995?)
(1, 0.2819716800611949?, 0.5310100564595692?)
(1, 1.883203505913526?, 0.5310100564595692?)
(1, 0.5310100564595692?, 1.883203505913526?)
(1, 0.5310100564595692?, 0.2819716800611949?)
"""
# The limit roots below are not in an exact ring for speed reasons.
hyper_lm = Set([])
for element in self.hyperbolic_elements(length):
EigVec = self._algebraic_eigenvectors_element(element) # THIS SHOULD BE FAST
new_lr = []
for ev in EigVec: # Making sure it is not a unimodular eigenvalue
if abs(ev[0].norm())-1 > 0.001: # The eigenvector with reverse eigenvalue is obtained with the inverse element
new_lr += [ev[1][0]]
for vect in new_lr:
if vect < 0: # Putting the eigenvector in the positive cone
vect = -vect
real_vect = vector([v.real() for v in vect])
real_vect.set_immutable()
if real_vect not in hyper_lm:
hyper_lm = hyper_lm.union(Set([real_vect]))
self._computed_limit_roots = self._computed_limit_roots.union(hyper_lm)
return hyper_lm
[docs] @cached_method
def limit_roots(self, length):
r"""
For Lorentzian Coxeter groups, return the limit roots coming from the
eigenvectors of infinite order elements of length ``length``
INPUT:
- ``length`` -- a non-negative integer
OUTPUT:
A set of vectors
EXAMPLES::
sage: from brocoli import *
sage: M1 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,oo,oo],[oo,oo,1,oo],[oo,oo,oo,1]])
sage: M2 = CoxeterMatrix([[1,3,3,3],[3,1,3,3],[3,3,1,3],[3,3,3,1]])
sage: M3 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,3,3],[oo,3,1,3],[oo,3,3,1]])
sage: M4 = CoxeterMatrix([[1,5,oo,oo],[5,1,3,3],[oo,3,1,3],[oo,3,3,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: GR2 = GeometricRepresentationCoxeterGroup(M2)
sage: GR3 = GeometricRepresentationCoxeterGroup(M3)
sage: GR4 = GeometricRepresentationCoxeterGroup(M4)
sage: [len(GR1.limit_roots(i)) for i in range(4)]
[0, 0, 6, 24]
sage: [len(GR2.limit_roots(i)) for i in range(4)]
[0, 0, 0, 4]
sage: [len(GR3.limit_roots(i)) for i in range(4)]
[0, 0, 3, 19]
sage: [len(GR4.limit_roots(i)) for i in range(4)] # long time (8 seconds)
[0, 0, 2, 19]
TESTS::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,4,oo,oo],[4,1,oo,oo],[oo,oo,1,-2],[oo,oo,-2,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: l_roots = GR.limit_roots(4) # long time (60 seconds)
There was an approximation problem with the element 1241
There was an approximation problem with the element 1421
There was an approximation problem with the element 1321
There was an approximation problem with the element 1231
"""
return self.parabolic_limit_roots(length).union(self.hyperbolic_limit_roots(length))
[docs] def visualize_limit_roots(self, list_lengths, list_orbits=[0], limit_type=2, size=4, color=(1, 0, 0), color_simplex=(0, 1, 0), isotropic=False, iso_color=(1, 1, 0)):
r"""
Return a graphical representation of limit roots
INPUT:
- ``list_lengths`` -- list of non-negative integers; they give the length of the infinite
order elements to be used to obtain limit roots
- ``list_orbits`` -- a list of non-negative integers (default: ``[0]``); they give the length
of the elements that will act on the limit roots to obtain even more of them
- ``limit_type`` -- an integer (default: ``2``);
- 2: Compute all limit roots
- 1: Compute hyperbolic limit roots
- anything else: Compute parabolic limit roots
- ``size`` -- a non-negative integer (default: ``4``); give the size of the point to draw
- ``color`` -- a color to give the limit roots (default: ``(1,0,0)``, red);
- ``color_simplex`` -- a color to give the limit roots (default: ``(0,1,0)``, green);
- ``isotropic`` -- a boolean (default: ``False``); whether to plot the
isotropic cone in the infinite case
- ``iso_color`` -- a color to give the isotropic cone (default: ``(1, 1, 0)``, yellow)
OUTPUT:
A Graphic Object
EXAMPLES::
sage: from brocoli import *
sage: DiAf = CoxeterMatrix([[1,oo],[oo,1]])
sage: GRDiAf = GeometricRepresentationCoxeterGroup(DiAf)
sage: img = GRDiAf.visualize_limit_roots([2]);img
Graphics object consisting of 2 graphics primitives
sage: DiHy = CoxeterMatrix([[1,-2],[-2,1]])
sage: GRDiHy = GeometricRepresentationCoxeterGroup(DiHy)
sage: img = GRDiHy.visualize_limit_roots([2]);img
Graphics object consisting of 3 graphics primitives
sage: M1 = CoxeterMatrix([[1,4,4],[4,1,4],[4,4,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: img = GR1.visualize_limit_roots([0,1,2]);img
Graphics object consisting of 3 graphics primitives
sage: img = GR1.visualize_limit_roots([3,4]);img # long time
Graphics object consisting of 21 graphics primitives
sage: img = GR1.visualize_limit_roots([3,4],[0,1,2]);img # long time
Graphics object consisting of 183 graphics primitives
sage: M2 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,oo,oo],[oo,oo,1,oo],[oo,oo,oo,1]])
sage: GR2 = GeometricRepresentationCoxeterGroup(M2)
sage: img = GR2.visualize_limit_roots([2,3]);img # long time
Graphics3d Object
sage: img = GR2.visualize_limit_roots([2,3],limit_type=1);img # long time
Graphics3d Object
sage: img = GR2.visualize_limit_roots([2,3],limit_type=0);img
Graphics3d Object
"""
if self._rank <= 4:
img = Graphics()
if not self.is_finite() and isotropic:
img += self.visualize_isotropic_cone(color=iso_color)
img += plot_simplex(self._rank, color_simplex)
for length in list_lengths:
for orbit in list_orbits:
img += self._visualize_limit_roots(length, orbit, limit_type, size, color)
return img
else:
raise NotImplementedError("visualization of rank >=5")
@cached_method
def _visualize_limit_roots(self, length, orbit=0, limit_type=2, size=4, color=(1, 0, 0)):
r"""
Return a Graphics object containing the limit roots coming from element
of length ``length`` acted upon by elements of length ``orbit``.
INPUT:
- ``length`` -- a non-negative integer
- ``orbit`` -- a non-negative integer (default: ``0``)
- ``limit_type`` -- an integer (default: ``2``)
- 2: Compute all limit roots
- 1: Compute hyperbolic limit roots
- anything else: Compute parabolic limit roots
- ``size`` -- a non-negative integer (default: ``4``); give the size of the point to draw
- ``color`` -- a color to give the limit roots (default: ``(1,0,0)``, red)
OUTPUT:
A Graphics object
EXAMPLES::
sage: from brocoli import *
sage: DiAf = CoxeterMatrix([[1,oo],[oo,1]])
sage: GRDiAf = GeometricRepresentationCoxeterGroup(DiAf)
sage: img = GRDiAf._visualize_limit_roots(2);img
Graphics object consisting of 1 graphics primitive
sage: DiHy = CoxeterMatrix([[1,-2],[-2,1]])
sage: GRDiHy = GeometricRepresentationCoxeterGroup(DiHy)
sage: img = GRDiHy._visualize_limit_roots(2);img
Graphics object consisting of 2 graphics primitives
sage: M1 = CoxeterMatrix([[1,4,4],[4,1,4],[4,4,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: img = GR1._visualize_limit_roots(3);img
Graphics object consisting of 6 graphics primitives
sage: img = GR1._visualize_limit_roots(4);img # long time (5 seconds)
Graphics object consisting of 12 graphics primitives
sage: M1 = CoxeterMatrix([[1,4,4],[4,1,4],[4,4,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: img = GR1.visualize_limit_roots([2]);img
Graphics object consisting of 3 graphics primitives
sage: img = GR1.visualize_limit_roots([3]);img
Graphics object consisting of 9 graphics primitives
sage: img = GR1.visualize_limit_roots([4]);img
Graphics object consisting of 15 graphics primitives
sage: img = GR1.visualize_limit_roots([3,4]);img
Graphics object consisting of 21 graphics primitives
"""
if self._rank == 1:
raise ValueError("The Coxeter group does not have limit roots")
elif self._rank <= 4:
img = Graphics()
if self._rank == 2:
projection_space = VectorSpace(RDF, self._rank)
else:
projection_space = VectorSpace(RDF, self._rank-1)
affine_basis = regular_simplex_vertices(self._rank - 1)
set_limit_roots = self._compute_orbit_limit_roots(length, orbit, limit_type)
for lr in set_limit_roots:
projected_lr = affinely_project_vector(lr, projection_space, affine_basis)
img += point(projected_lr, size=size, rgbcolor=color)
if self._rank != 4:
img.set_aspect_ratio(1)
img.SHOW_OPTIONS['axes'] = False
return img
else:
raise NotImplementedError("visualization of rank >=5 is not available")
@cached_method
def _compute_orbit_limit_roots(self, base_length, orbit_length, limit_type=2):
r"""
Return the image of the action of element of length ``orbit_length`` on
the limit roots coming from elements of length ``base_length``
INPUT:
- ``base_length`` -- a non-negative integer
- ``orbit_length`` -- a non-negative integer
- ``limit_type`` -- an integer (default: ``2``)
- 2: Compute all limit roots
- 1: Compute hyperbolic limit roots
- anything else: Compute parabolic limit roots
OUTPUT:
A set of limit roots
EXAMPLES::
sage: from brocoli import *
sage: M1 = CoxeterMatrix([[1,4,4],[4,1,4],[4,4,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: GR1._compute_orbit_limit_roots(3,0)[0]
(1, 3.546455444684995?, 1.883203505913526?)
sage: GR1._compute_orbit_limit_roots(3,1)[-1]
(1.0, 0.5310100564595692, 0.2819716800611949)
sage: GR1._compute_orbit_limit_roots(3,2)[-1]
(6.678697326996895, 3.546455444684995, 1.8832035059135257)
sage: M2 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,3,3],[oo,3,1,3],[oo,3,3,1]])
sage: GR2 = GeometricRepresentationCoxeterGroup(M2)
sage: len(GR2._compute_orbit_limit_roots(2,0))
3
sage: len(GR2._compute_orbit_limit_roots(2,0,0))
3
sage: len(GR2._compute_orbit_limit_roots(2,0,1))
0
sage: len(GR2._compute_orbit_limit_roots(3,0))
19
sage: GR2._compute_orbit_limit_roots(3,0,0)
{(0.0, 1.0, 1.0, 1.0)}
sage: len(GR2._compute_orbit_limit_roots(3,0,1))
18
sage: len(GR2._compute_orbit_limit_roots(3,1))
74
sage: GR2._compute_orbit_limit_roots(3,1,0)
{(6.0, 1.0, 1.0, 1.0), (0.0, 1.0, 1.0, 1.0)}
sage: len(GR2._compute_orbit_limit_roots(3,1,1))
72
TODO:
Solve the approximation issue
"""
if limit_type == 2:
limit_roots = self.limit_roots(base_length)
elif limit_type == 1:
limit_roots = self.hyperbolic_limit_roots(base_length)
else:
limit_roots = self.parabolic_limit_roots(base_length)
orbit_lm = Set([])
if orbit_length == 0:
return limit_roots
else:
for lm in limit_roots:
for element in self.elements(orbit_length):
new_lm = self._real_element(element)*lm
if new_lm < 0:
new_lm = -new_lm
approx_problem = False
index = 0
while not approx_problem and index < self._rank:
if new_lm[index] < 0: # SOLVE THIS ISSUE HERE
approx_problem = True
# print "There was an approximation problem with element:", self._matrix_to_word[element], "and limit root:", lm
index += 1
if not approx_problem:
new_lm.set_immutable()
if new_lm not in orbit_lm:
orbit_lm = orbit_lm.union(Set([new_lm]))
self._computed_limit_roots = self._computed_limit_roots.union(orbit_lm)
return orbit_lm
[docs] def visualize_roots(self, list_depths, size=4, color=(1, 0, 0), color_simplex=(0, 1, 0), isotropic=True, iso_color=(1, 1, 0)):
r"""
Return a graphical representation of roots of depths in ``list_depths``
INPUT:
- ``list_depths`` -- a list of non-negative integers; give the depths of
the roots to plot
- ``size`` -- a non-negative integer (default: ``4``); give the size of the point to draw
- ``color`` -- a color to give the roots (default: ``(1,0,0)``, red)
- ``color_simplex`` -- a color to give the simplex (default: ``(0,1,0)``, green)
- ``isotropic`` -- a boolean (default: ``True``); whether to plot the
isotropic cone in the infinite case
- ``iso_color`` -- a color to give the isotropic cone (default: ``(1, 1, 0)``, yellow)
OUTPUT:
A Graphic Object
EXAMPLES:
Rank 2 examples::
sage: from brocoli import *
sage: A2 = CoxeterMatrix([[1,3],[3,1]])
sage: GRA2 = GeometricRepresentationCoxeterGroup(A2)
sage: img = GRA2.visualize_roots(range(2));img
Graphics object consisting of 4 graphics primitives
sage: B2 = CoxeterMatrix([[1,4],[4,1]])
sage: GRB2 = GeometricRepresentationCoxeterGroup(B2)
sage: img = GRB2.visualize_roots(range(2));img
Graphics object consisting of 5 graphics primitives
sage: DiAf = CoxeterMatrix([[1,oo],[oo,1]])
sage: GRDiAf = GeometricRepresentationCoxeterGroup(DiAf)
sage: img = GRDiAf.visualize_roots(range(5));img
Graphics object consisting of 12 graphics primitives
sage: DiHy = CoxeterMatrix([[1,-2],[-2,1]])
sage: GRDiHy = GeometricRepresentationCoxeterGroup(DiHy)
sage: img = GRDiHy.visualize_roots(range(2));img
Graphics object consisting of 8 graphics primitives
Rank 3 examples::
sage: A3 = CoxeterMatrix([[1,3,2],[3,1,3],[2,3,1]])
sage: GRA3 = GeometricRepresentationCoxeterGroup(A3)
sage: img = GRA3.visualize_roots(range(3));img
Graphics object consisting of 9 graphics primitives
sage: B3 = CoxeterMatrix([[1,4,2],[4,1,3],[2,3,1]])
sage: GRB3 = GeometricRepresentationCoxeterGroup(B3)
sage: img = GRB3.visualize_roots(range(4));img
Graphics object consisting of 12 graphics primitives
sage: H3 = CoxeterMatrix([[1,5,2],[5,1,3],[2,3,1]])
sage: GRH3 = GeometricRepresentationCoxeterGroup(H3)
sage: img = GRH3.visualize_roots(range(7));img
Graphics object consisting of 18 graphics primitives
sage: J3 = CoxeterMatrix([[1,6,2],[6,1,3],[2,3,1]])
sage: GRJ3 = GeometricRepresentationCoxeterGroup(J3)
sage: img = GRJ3.visualize_roots(range(10));img
Graphics object consisting of 40 graphics primitives
sage: K3 = CoxeterMatrix([[1,7,2],[7,1,3],[2,3,1]])
sage: GRK3 = GeometricRepresentationCoxeterGroup(K3)
sage: img = GRK3.visualize_roots(range(10));img
Graphics object consisting of 55 graphics primitives
Rank 4 examples::
sage: A4 = CoxeterMatrix([[1,3,2,2],[3,1,3,2],[2,3,1,3],[2,2,3,1]])
sage: GRA4 = GeometricRepresentationCoxeterGroup(A4)
sage: img = GRA4.visualize_roots(range(4));img
Graphics3d Object
sage: C4t = CoxeterMatrix([[1,4,2,2],[4,1,3,2],[2,3,1,4],[2,2,4,1]])
sage: GRC4t = GeometricRepresentationCoxeterGroup(C4t)
sage: img = GRC4t.visualize_roots(range(6));img # long time
Graphics3d Object
sage: M = CoxeterMatrix([[1,oo,oo,oo],[oo,1,3,3],[oo,3,1,3],[oo,3,3,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: img = GR.visualize_roots(range(4));img # long time
Graphics3d Object
"""
img = Graphics()
img += plot_simplex(self._rank, color_simplex)
if not self.is_finite() and isotropic:
img += self.visualize_isotropic_cone(color=iso_color)
if self._rank <= 4:
for depth in list_depths:
img += self._visualize_roots_depth(depth, size, color)
if self._rank <= 3:
img.set_aspect_ratio(1)
img.SHOW_OPTIONS['axes'] = False
return img
else:
raise NotImplementedError("visualization of rank >=5 is not available")
@cached_method
def _visualize_roots_depth(self, depth, size=4, color=(1, 0, 0)):
r"""
Return a graphical representation of roots of depth ``depth``
INPUT:
- ``depth`` -- a non-negative integer
- ``size`` -- a non-negative integer (default: ``4``); give the size of the point to draw
- ``color`` -- a color to give the points (default: ``(1,0,0)``, red)
OUTPUT:
A Graphic Object
EXAMPLES::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,oo,oo,oo],[oo,1,3,3],[oo,3,1,3],[oo,3,3,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: GR._visualize_roots_depth(0)
Graphics3d Object
sage: GR._visualize_roots_depth(1)
Graphics3d Object
"""
if self._rank == 2:
real_roots = [self._rdf_vector(root) for root in self.roots(depth)]
img = Graphics()
for rr in real_roots:
img += arrow((0, 0), tuple(rr), color=color)
return img
elif self._rank <= 4:
projection_space = VectorSpace(RDF, self._rank-1)
affine_basis = regular_simplex_vertices(self._rank - 1)
img = Graphics()
for r in self.roots(depth):
img += point(list(affinely_project_vector(r, projection_space, affine_basis)), size=size, rgbcolor=color)
return img
else:
raise NotImplementedError("visualization of rank >=5 roots.")
[docs] def visualize_weights(self, list_depths, space=False, size=4, color=(0, 0, 1), color_simplex=(0, 1, 0), isotropic=True, iso_color=(1, 1, 0)):
r"""
Return a graphical representation of weights. Some weights may not have
an affine representation and will not be drawn.
INPUT:
- ``list_depths`` -- a list of non-negative integers; give the depths of
the weights to plot
- ``space`` -- a boolean (default: ``False``); whether to plot only
space weights
- ``size`` -- a non-negative integer (default: ``4``); give the size of the point to draw
- ``color`` -- a color to give the weights (default: ``(1,0,0)``, blue)
- ``color_simplex`` -- a color to give the limit roots (default: ``(0,1,0)``, green)
- ``isotropic`` -- a boolean (default: ``True``); whether to plot the
isotropic cone in the infinite case
- ``iso_color`` -- a color to give the isotropic cone (default: ``(1, 1, 0)``, yellow)
OUTPUT:
A Graphic Object
EXAMPLES:
Rank 2 examples::
sage: from brocoli import *
sage: A2 = CoxeterMatrix([[1,3],[3,1]])
sage: GRA2 = GeometricRepresentationCoxeterGroup(A2)
sage: GRA2.visualize_weights(range(3))
Graphics object consisting of 7 graphics primitives
sage: DiHy = CoxeterMatrix([[1,-2],[-2,1]])
sage: GRDiHy = GeometricRepresentationCoxeterGroup(DiHy)
sage: GRDiHy.visualize_weights(range(3))
Graphics object consisting of 10 graphics primitives
Rank 3 examples::
sage: A3 = CoxeterMatrix([[1,3,2],[3,1,3],[2,3,1]])
sage: GRA3 = GeometricRepresentationCoxeterGroup(A3)
sage: GRA3.visualize_weights(range(5))
Graphics object consisting of 15 graphics primitives
sage: M = CoxeterMatrix([[1,-10/9,oo],[-10/9,1,4],[oo,4,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: GR.visualize_weights(range(4))
Graphics object consisting of 31 graphics primitives
Rank 4 examples::
sage: A4 = CoxeterMatrix([[1,3,2,2],[3,1,3,2],[2,3,1,3],[2,2,3,1]])
sage: GRA4 = GeometricRepresentationCoxeterGroup(A4)
sage: GRA4.visualize_weights(range(7)) # Not all weights are drawn
Graphics3d Object
sage: C4t = CoxeterMatrix([[1,4,2,2],[4,1,3,2],[2,3,1,4],[2,2,4,1]])
sage: GRC4t = GeometricRepresentationCoxeterGroup(C4t)
sage: GRC4t.visualize_weights(range(5))
Traceback (most recent call last):
...
ValueError: The weights can not be visualized affinely; the bilinear form is degenerate
sage: QF.<a> = QuadraticField(3/2)
sage: Cyl = CoxeterMatrix([[1,-a,2,2],[-a,1,3,2],[2,3,1,-a],[2,2,-a,1]])
sage: GRCyl = GeometricRepresentationCoxeterGroup(Cyl)
sage: GRCyl.visualize_weights(range(5))
Traceback (most recent call last):
...
ValueError: The weights can not be visualized affinely; the bilinear form is degenerate
sage: M = CoxeterMatrix([[1,oo,oo,oo],[oo,1,3,3],[oo,3,1,3],[oo,3,3,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: GR.visualize_weights(range(3)) # long time
Graphics3d Object
sage: GR.visualize_weights(range(3),True) # This one has less weights
Graphics3d Object
WARNING:
In certain cases, some weights have a sum of coordinates equal to
zero, they are ignored by the plotting function.
"""
if self.signature()[2] != 0:
raise ValueError("The weights can not be visualized affinely; the bilinear form is degenerate")
img = Graphics()
img += plot_simplex(self._rank, color_simplex)
if not self.is_finite() and isotropic:
img += self.visualize_isotropic_cone(color=iso_color)
if self._rank <= 4:
for depth in list_depths:
img += self._visualize_weights_depth(depth, space, size, color)
if self._rank <= 3:
img.set_aspect_ratio(1)
img.SHOW_OPTIONS['axes'] = False
return img
else:
raise NotImplementedError("visualization of rank >=5 is not available")
@cached_method
def _visualize_weights_depth(self, depth, space=False, size=4, color=(0, 0, 1)):
r"""
Return a graphical representation of weights of depth ``depth``
INPUT:
- ``depth`` -- a non-negative integer
- ``space`` -- a boolean (default: ``False``); whether to plot only space
weights
- ``size`` -- a non-negative integer (default: ``4``); give the size of the point to draw
- ``color`` -- a color to give the weights (default: ``(0,0,1)``, blue)
OUTPUT:
A Graphic Object
EXAMPLES::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,oo,oo,oo],[oo,1,3,3],[oo,3,1,3],[oo,3,3,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: GR._visualize_weights_depth(0)
Graphics3d Object
sage: GR._visualize_weights_depth(1)
Graphics3d Object
"""
if self._rank == 2:
if space:
real_weights = [self._rdf_vector(weight) for weight in self.space_weights(depth)]
else:
real_weights = [self._rdf_vector(weight) for weight in self.weights(depth)]
img = Graphics()
for rr in real_weights:
img += arrow((0, 0), tuple(rr), color=color)
return img
elif self._rank <= 4:
projection_space = VectorSpace(RDF, self._rank-1)
affine_basis = regular_simplex_vertices(self._rank - 1)
img = Graphics()
if space:
plotted_weights = self.space_weights(depth)
else:
plotted_weights = self.weights(depth)
for w in plotted_weights:
try:
img += point(list(affinely_project_vector(w, projection_space, affine_basis)), size=size, rgbcolor=color)
except ValueError: # Not every weight can be affinely projected!
pass
return img
else:
raise NotImplementedError("visualization of rank >=5 roots.")
[docs] @cached_method
def visualize_tits_cone(self, depth, frame_color=(1, 0, 0), face_color=(0, 0, 1), color_simplex=(0, 1, 0), isotropic=False, iso_color=(1, 1, 0)):
r"""
Return a graphical representation of the Tits cone up to weights of
depth ``depth``.
INPUT:
- ``depth`` -- a non-negative integer
- ``frame_color`` -- a color (default= ``(1,0,0)`` red); color of the
edges
- ``face_color`` -- a color (default= ``(0, 0, 1)``); color of the faces
- ``color_simplex`` -- a color to give the limit roots (default: ``(0,1,0)``, green)
- ``isotropic`` -- a boolean (default: ``False``); whether to plot the isotropic cone
- ``iso_color`` -- a color to give the isotropic cone (default: ``(1, 1, 0)``, yellow)
OUTPUT:
A Graphic Object
EXAMPLES::
sage: from brocoli import *
sage: M1 = CoxeterMatrix([[1,4,4],[4,1,4],[4,4,1]])
sage: GR1 = GeometricRepresentationCoxeterGroup(M1)
sage: GR1.visualize_tits_cone(2)
Graphics object consisting of 13 graphics primitives
sage: M2 = CoxeterMatrix([[1,-2,-2],[-2,1,-2],[-2,-2,1]])
sage: GR2 = GeometricRepresentationCoxeterGroup(M2)
sage: GR2.visualize_tits_cone(2)
Graphics object consisting of 16 graphics primitives
sage: M3 = CoxeterMatrix([[1,oo,oo,oo],[oo,1,3,3],[oo,3,1,3],[oo,3,3,1]])
sage: GR3 = GeometricRepresentationCoxeterGroup(M3)
sage: GR3.visualize_tits_cone(3) # long time (5 seconds)
Graphics3d Object
"""
if self.is_finite():
raise ValueError("Can not draw the Tits cone for finite Coxeter groups")
elif self.signature()[2] != 0:
raise ValueError("The Tits cone cannot be visualized affinely; the bilinear form is degenerate")
elif self._rank <= 4:
projection_space = VectorSpace(RDF, self._rank - 1)
affine_basis = regular_simplex_vertices(self._rank - 1)
proj_weights = []
for dp in range(depth + 1):
for weight in self.weights(dp):
proj_weights += [affinely_project_vector(weight, projection_space, affine_basis)]
P = Polyhedron(vertices=proj_weights, base_ring=RDF)
img = plot_simplex(self._rank, color_simplex)
if isotropic:
img += self.visualize_isotropic_cone(color=iso_color)
img += P.render_wireframe(color=frame_color)
img += P.render_solid(alpha=0.5, color=face_color)
return img
else:
raise NotImplementedError("rank >=5")
[docs] @cached_method
def visualize_isotropic_cone(self, color=(1, 1, 0), color_simplex=(0, 1, 0), simplex=True):
r"""
Return a graphical representation of the isotropic cone
INPUT:
- ``color`` -- a color (default= ``(1,1,0)``)
- ``color_simplex`` -- a color to give the limit roots (default: ``(0,1,0)``, green)
- ``simplex`` -- a boolean (default: ``True``); whether to plot the affine simplex
OUTPUT:
A Graphic Object
EXAMPLES::
sage: from brocoli import *
sage: A3 = CoxeterMatrix([[1,3,2],[3,1,3],[2,3,1]])
sage: GRA3 = GeometricRepresentationCoxeterGroup(A3)
sage: GRA3.visualize_isotropic_cone()
Traceback (most recent call last):
...
ValueError: The bilinear form is positive definite
sage: C4t = CoxeterMatrix([[1,4,2,2],[4,1,3,2],[2,3,1,4],[2,2,4,1]])
sage: GRC4t = GeometricRepresentationCoxeterGroup(C4t)
sage: GRC4t.visualize_isotropic_cone()
Graphics3d Object
sage: M = CoxeterMatrix([[1,oo,oo,oo],[oo,1,3,3],[oo,3,1,3],[oo,3,3,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: GR.visualize_isotropic_cone()
Graphics3d Object
"""
if self.is_finite():
raise ValueError("The bilinear form is positive definite")
elif self.is_affine():
pt = self.bilinear_form().kernel().basis()[0]
rdf_pt = self._rdf_vector(pt) # Done to avoid problem with norm in CyclotomicField
pt_norm = 5 * rdf_pt/rdf_pt.norm()
if self._rank == 2:
image = arrow((0, 0), list(pt_norm), color=color)
return image
elif self._rank <= 4:
projection_space = VectorSpace(RDF, self._rank-1)
affine_basis = regular_simplex_vertices(self._rank - 1)
center = affinely_project_vector(pt_norm, projection_space, affine_basis)
image = point(list(center), size=8, rgbcolor=color)
if self._rank == 3:
image += circle(list(center), 0.01, rgbcolor=(0, 0, 0))
image += plot_simplex(self._rank, color_simplex)
return image
else:
raise NotImplementedError("rank >=5")
elif self._rank <= 4:
return self._draw_isotropic_cone_surface(color) + plot_simplex(self._rank, color_simplex)
else:
raise NotImplementedError("rank >=5")
def _draw_isotropic_cone_surface(self, color=(1, 1, 0)):
r"""
Return a graphical object representing the isotropic cone for infinite
non-affine Coxeter groups of rank 3 or 4
INPUT:
- ``color`` -- a color (default= ``(1,1,0)``)
OUTPUT:
A Graphic Object
EXAMPLES::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,oo,oo,oo],[oo,1,3,3],[oo,3,1,3],[oo,3,3,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: GR._draw_isotropic_cone_surface()
Graphics3d Object
"""
if self._rank == 4:
eq_lc = self._equation_light_cone_rank4()
image = implicit_plot3d(eq_lc == 0, (a, -0.5, 2.1), (b, -0.5, 2.1), (c, -0.5, 2.1), aspect_ratio=1, rgbcolor=color, opacity=0.2)
return image
elif self._rank == 3:
eq_lc = self._equation_light_cone_rank3()
image = implicit_plot(eq_lc == 0, (a, -1, 2), (b, -1, 2), aspect_ratio=1, color=color)
image.set_aspect_ratio(1)
image.SHOW_OPTIONS['axes'] = False
return image
elif self._rank == 2:
eq_lc = self._equation_light_cone_rank2()
values = [soln.values()[0] for soln in solve(eq_lc == 0, a, solution_dict=True)]
v1 = 5 * vector(RDF, [1-values[0], values[0]])
v2 = 5 * vector(RDF, [1-values[1], values[1]])
image = arrow((0, 0), tuple(v1), color=color)
image += arrow((0, 0), tuple(v2), color=color)
return image
else:
raise ValueError("can not draw the isotropic cone for rank 1 or rank >=5")
def _equation_light_cone_rank4(self):
r"""
Return the implicit equation of the isotropic cone in the affine
3-space spanned by the simple roots
OUTPUT:
A symbolic expression
EXAMPLES::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,oo,oo,oo],[oo,1,3,3],[oo,3,1,3],[oo,3,3,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: GR._equation_light_cone_rank4().coefficients() # abs tol 1e-15
[[0.08333333333333338*sqrt(2)*b*c + 0.8333333333333333*b^2 - 0.33333333333333337*sqrt(3)*sqrt(2)*c + 0.7916666666666666*c^2 - 0.6666666666666666*sqrt(3)*b + 1.0,
0],
[0.08333333333333331*sqrt(3)*sqrt(2)*c + 0.16666666666666663*sqrt(3)*b - 2.0,
1],
[1.0, 2]]
"""
a, b, c = var('a b c')
v = vector([1-a-b-c, a, b, c])
vect = v*self.bilinear_form(False)*v.column()
my_eq = vect[0]
eq_mod = my_eq.subs(b=b/sqrt(3))
eq_mod = eq_mod.subs(a=a/2)
eq_mod = eq_mod.subs(c=c*(sqrt(3)/(2*sqrt(2))))
eq_mod = eq_mod.subs(a=a-b/sqrt(3))
eq_mod = eq_mod.subs(b=b-c/(2*sqrt(2)))
eq_mod = eq_mod.subs(a=a-c*(sqrt(3)/(2*sqrt(2))))
return eq_mod
def _equation_light_cone_rank3(self):
r"""
Return the implicit equation of the isotropic cone in the affine plane
spanned by the simple roots
OUTPUT:
A symbolic expression
EXAMPLES::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,-10/9,oo],[-10/9,1,4],[oo,4,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: GR._equation_light_cone_rank3().coefficients() # abs tol 1e-15
[[0.8838504085436638*b^2 - 0.6296296296296295*sqrt(3)*b + 1.0, 0],
[0.09763107293781749*sqrt(3)*b - 2.111111111111111, 1],
[1.0555555555555556, 2]]
"""
a, b = var('a b')
v = vector([1-a-b, a, b])
vect = v*self.bilinear_form(False)*v.column()
my_eq = vect[0]
eq_mod = my_eq.subs(b=b/sqrt(3))
eq_mod = eq_mod.subs(a=a/2)
eq_mod = eq_mod.subs(a=a-b/sqrt(3))
return eq_mod
def _equation_light_cone_rank2(self):
r"""
Return the implicit equation of the isotropic cone in the affine plane
spanned by the simple roots
OUTPUT:
A symbolic expression
EXAMPLES::
sage: from brocoli import *
sage: DiHy = CoxeterMatrix([[1,-2],[-2,1]])
sage: GRDiHy = GeometricRepresentationCoxeterGroup(DiHy)
sage: GRDiHy._equation_light_cone_rank2()
(3.0*a - 2.0)*a - (a - 1)*(-3.0*a + 1.0)
"""
a = var('a')
v = vector([1-a, a])
vect = v*self.bilinear_form(False)*v.column()
my_eq = vect[0]
return my_eq
[docs] def tikz_picture_rank3(self, range_roots, range_limitroots, range_orbit,
range_weights, range_limitdir, range_hyperplane, limit_type=2, isotropic=True):
r"""
This function returns a LaTeX expression containing a tikzpicture of a rank
3 root system.
INPUT:
- ``range_roots`` -- list of non-negative integers; the depths of roots to print
- ``range_limitroots`` -- list of non-negative integers; give the limit roots associated to elements
with lengths in `range_limitroots`
- ``range_orbit`` -- list of non-negative integers; lengths by which we act on roots and limit roots to obtain more
(be careful not to put to many)
- ``range_weights`` -- list of non-negative integers; depth of weights to print
- ``range_limitdir`` -- list of non-negative integers; lengths of elements of which we print the
possible associated limit direction (eigenspace of dimension 1)
- ``range_hyperplane`` -- list of non-negative integers; depths of roots for which we print the
corresponding hyperplane inside the isotropic cone
- ``limit_type`` -- 0, 1, 2 (default: ``2``); ``0`` : Parabolic, ``1`` : Hyperbolic, or ``2`` : All
- ``isotropic`` -- Boolean (default: ``True``); whether to draw the isotropic cone
OUTPUT:
A string containing a tikzpicture.
EXAMPLES::
sage: from brocoli import *
sage: M = CoxeterMatrix([[1,-10/9,oo],[-10/9,1,4],[oo,4,1]])
sage: GR = GeometricRepresentationCoxeterGroup(M)
sage: tikz_pic = GR.tikz_picture_rank3(range(3),[2,3],[0,1],range(3),[2,3],range(2)) # long time (8 seconds)
sage: print("\n".join(tikz_pic.splitlines()[:10])) # long time
\begin{tikzpicture}
[scale=1,
isotropcone/.style={red,line join=round,thick},
root/.style={blue},
simpleroot/.style={black},
dihedral_roots/.style={blue},
limit/.style={fill=red,draw=black,diamond},
limdir/.style={fill=orange,draw=black,diamond},
weight/.style={fill=green,draw=black,diamond},
rotate=0]
sage: print("\n".join(tikz_pic.splitlines()[69:74])) # long time
\node[limit,inner sep=1pt] at (2.79384317461,0.308476464594) {};
\node[limit,inner sep=1pt] at (1.6411385843,2.19814939939) {};
\node[limit,inner sep=1pt] at (1.11895251609,0.356038263561) {};
\node[limit,inner sep=1pt] at (2.35114381286,2.06826137522) {};
\node[limit,inner sep=1pt] at (0.900982021125,0.792578618065) {};
"""
if self._rank != 3:
raise ValueError("The representation does not have rank 3")
tikz_pic = '\\begin{tikzpicture}\n'
tikz_pic += '\t[scale=1,\n'
tikz_pic += '\t isotropcone/.style={red,line join=round,thick},\n'
tikz_pic += '\t root/.style={blue},\n'
tikz_pic += '\t simpleroot/.style={black},\n'
tikz_pic += '\t dihedral_roots/.style={blue},\n'
tikz_pic += '\t limit/.style={fill=red,draw=black,diamond},\n'
tikz_pic += '\t limdir/.style={fill=orange,draw=black,diamond},\n'
tikz_pic += '\t weight/.style={fill=green,draw=black,diamond},\n'
tikz_pic += '\t rotate=0]\n\n\n'
tikz_pic += '\def\size{0.0125}\n'
tikz_pic += '\def\simplesize{0.025}\n\n'
projection_space = VectorSpace(RDF, 2)
affine_basis = regular_simplex_vertices(2)
if self.signature() == [2, 0, 1] and isotropic:
bf = self.bilinear_form()
pt = bf.kernel().basis()[0]
center = affinely_project_vector(vector(pt), projection_space, affine_basis)
tikz_pic += '\\fill[red] (%s,%s) circle (\simplesize);\n\n' % (2*center[0], 2*center[1])
elif isotropic:
from sage.plot.misc import setup_for_eval_on_grid
eq_lc = self._equation_light_cone_rank3()
g, ranges = setup_for_eval_on_grid([eq_lc], [[a, 0, 4], [b, 0, 4]], 200)
g = g[0]
pt_sequence = []
centerval = 0
center = None
for i in xsrange(-0.5, 4, step=0.005):
for j in xsrange(-0.5, 4, step=0.005):
val = g(i, j)
if abs(val) < 0.00125 and val < 0:
pt_sequence += [(i, j)]
if val < centerval:
centerval = val
center = (i, j)
segments = [(p[0] - center[0], p[1] - center[1]) for p in pt_sequence]
segments_pospos = filter(lambda x: x[0] > 0 and x[1] > 0, segments)
segments_negpos = filter(lambda x: x[0] < 0 and x[1] > 0, segments)
segments_negneg = filter(lambda x: x[0] < 0 and x[1] <= 0, segments)
segments_posneg = filter(lambda x: x[0] > 0 and x[1] <= 0, segments)
segpospos_tan = [(max(arctan(p[1]/p[0]), -arctan(p[1]/p[0])), p) for p in segments_pospos]
segnegpos_tan = [(min(arctan(p[1]/p[0]), -arctan(p[1]/p[0])), p) for p in segments_negpos]
segnegneg_tan = [(max(arctan(p[1]/p[0]), -arctan(p[1]/p[0])), p) for p in segments_negneg]
segposneg_tan = [(min(arctan(p[1]/p[0]), -arctan(p[1]/p[0])), p) for p in segments_posneg]
segpospos_tan.sort()
segnegpos_tan.sort()
segnegneg_tan.sort()
segposneg_tan.sort()
Sorted_pt = [(p[1][0] + center[0], p[1][1] + center[1]) for p in segpospos_tan] + \
[(p[1][0] + center[0], p[1][1] + center[1]) for p in segnegpos_tan] + \
[(p[1][0] + center[0], p[1][1] + center[1]) for p in segnegneg_tan] + \
[(p[1][0] + center[0], p[1][1] + center[1]) for p in segposneg_tan]
impl = '%% The isotropic curve:\n\n\\draw[q] '
for p in Sorted_pt:
x_coord = n(2*p[0])
y_coord = n(2*p[1])
if abs(x_coord) < 0.001:
x_coord = n(0)
if abs(y_coord) < 0.001:
y_coord = n(0)
impl += '(%s,%s) -- ' % (x_coord.n(digits=3), y_coord.n(digits=3))
impl += 'cycle;\n\n'
tikz_pic += impl
simple = ['\\alpha', '\\beta', '\\gamma']
name = ['a', 'b', 'g']
place = ['left', 'right', 'above']
si_root = 0
for i in range_roots:
tikz_pic += '%% Roots of depth= %i\n' % i
for j in self.roots(i):
projected_root = affinely_project_vector(self._rdf_vector(j), projection_space, affine_basis)
line_sgmt = ''
if i == 0:
line_sgmt += '\\node[label=%s :{$%s$}] (%s) at (%s,%s) {};\n\\fill[simpleroot] (%s,%s) circle (\simplesize);\n' % (place[si_root], simple[si_root], name[si_root], 2*projected_root[0], 2*projected_root[1], 2*projected_root[0], 2*projected_root[1])
si_root += 1
else:
if abs(sqrt(3)*projected_root[0]-projected_root[1]) < 0.01 or abs(-sqrt(3)*projected_root[0]+2*sqrt(3)-projected_root[1]) < 0.01 or abs(projected_root[1]) < 0.01:
line_sgmt += '\\fill[root] (%s,%s) circle (\simplesize);\n' % (max(2*projected_root[0], 0.01), max(2*projected_root[1], 0.01))
else:
line_sgmt += '\\fill[root] (%s,%s) circle (\size);\n' % (max(2*projected_root[0], 0.01), max(2*projected_root[1], 0.01))
tikz_pic += line_sgmt
tikz_pic += '\n'
if i == 0:
tikz_pic += "\draw[green!75!black] (%s) -- (%s) -- (%s) -- (%s);\n\n" % (name[0], name[1], name[2], name[0])
for i in range_limitroots:
for j in range_orbit:
tikz_pic += '%% Limit roots of of type= %i and length= %i and orbit= %i\n\n' % (limit_type, i, j)
limit_roots = self._compute_orbit_limit_roots(i, j, limit_type)
for lm in limit_roots:
projected_lm = affinely_project_vector(lm, projection_space, affine_basis)
line_sgmt = ''
if abs(sqrt(3)*projected_lm[0]-projected_lm[1]) < 0.01 or abs(-sqrt(3)*projected_lm[0]+2*sqrt(3) - projected_lm[1]) < 0.01 or abs(projected_lm[1]) < 0.01:
line_sgmt += '\\node[limit,inner sep=1.5pt] at (%s,%s) {};\n' % (max(2*projected_lm[0], 0.01), max(2*projected_lm[1], 0.01))
else:
line_sgmt += '\\node[limit,inner sep=1pt] at (%s,%s) {};\n' % (max(2*projected_lm[0], 0.01), max(2*projected_lm[1], 0.01))
tikz_pic += line_sgmt
tikz_pic += '\n'
limdir_set = Set([])
x, y, z = var('x y z')
p = var('p')
var_vector = vector([x, y, z])
positive_cone = Polyhedron(rays=[[1, 0, 0], [0, 1, 0], [0, 0, 1]])
bilin_form = self.bilinear_form()
for depth in range_hyperplane:
tikz_pic += '%% Hyperplanes of depth= %i\n' % depth
for root in self.roots(depth):
equation1 = (root*bilin_form*var_vector.column())[0]
hyperplane = Polyhedron(ieqs=[[0, RDF(equation1.coefficient(x)), RDF(equation1.coefficient(y)), RDF(equation1.coefficient(z))], [0, RDF(-equation1.coefficient(x)), RDF(-equation1.coefficient(y)), RDF(-equation1.coefficient(z))]], base_ring=RDF)
positive_hyper = positive_cone.intersection(hyperplane)
ray0, ray1 = positive_hyper.rays()
affine_span = (1 - p) * vector(ray0) + p * vector(ray1)
isotropic_eq = (affine_span * bilin_form * affine_span.column())[0]
isotropic_rays = [(1 - i) * vector(ray0) + i * vector(ray1) for i in [(isotropic_eq.solve(p)[0]).rhs(), (isotropic_eq.solve(p)[1]).rhs()]]
affine_isorays = [affinely_project_vector(pt, projection_space, affine_basis) for pt in isotropic_rays]
tikz_pic += '\\draw[thin] ($(%s,%s)!3!(%s,%s)$) -- ($(%s,%s)!3!(%s,%s)$);\n' % (2*affine_isorays[0][0].n(digits=6), 2*affine_isorays[0][1].n(digits=6), 2*affine_isorays[1][0].n(digits=6), 2*affine_isorays[1][1].n(digits=6), 2*affine_isorays[1][0].n(digits=6), 2*affine_isorays[1][1].n(digits=6), 2*affine_isorays[0][0].n(digits=6), 2*affine_isorays[0][1].n(digits=6))
tikz_pic += '\n'
tikz_pic += '\n'
for i in range_limitdir:
tikz_pic += '%% Limit directions of of type= %i and length= %i\n\n' % (limit_type, i)
hyp_elts = self.hyperbolic_elements(i)
for hyper_el in hyp_elts:
eig_vectors = self._algebraic_eigenvectors_element(hyper_el)
uni_mod_vect = [ev[1][0] for ev in eig_vectors if abs(ev[0]) == 1]
new_ld = uni_mod_vect[0]
new_ld.set_immutable()
if new_ld not in limdir_set:
limdir_set = limdir_set.union(Set([new_ld]))
try:
projected_ld = affinely_project_vector(new_ld, projection_space, affine_basis)
line_sgmt = '\\node[limdir,inner sep=1pt] at (%s,%s) {}; pourc %s\n' % (2*projected_ld[0], 2*projected_ld[1], str(self.matrix_to_word(hyper_el)))
tikz_pic += line_sgmt
except ValueError:
pass
tikz_pic += '\n'
for i in range_weights:
tikz_pic += '%% Weights of length= %i\n' % i
for j in self.weights(i):
proj_weight = affinely_project_vector(j, projection_space, affine_basis)
line_sgmt = ''
if abs(2*proj_weight[0]) < 0.01:
x_coord = 0.01
else:
x_coord = 2*proj_weight[0]
if abs(2*proj_weight[1]) < 0.01:
y_coord = 0.01
else:
y_coord = 2*proj_weight[1]
if abs(sqrt(3)*proj_weight[0]-proj_weight[1]) < 0.01 or abs(-sqrt(3)*proj_weight[0]+2*sqrt(3)-proj_weight[1]) < 0.01 or abs(proj_weight[1]) < 0.01:
line_sgmt += '\\node[weight,inner sep=1.5pt] at (%s,%s) {};\n' % (x_coord, y_coord)
else:
line_sgmt += '\\node[weight,inner sep=1pt] at (%s,%s) {};\n' % (x_coord, y_coord)
tikz_pic += line_sgmt
tikz_pic += '\n'
tikz_pic += '\n\n\\end{tikzpicture}\n\n'
return LatexExpr(tikz_pic)