ClusterSeed¶
A cluster seed is a pair \((B,\mathbf{x})\) with \(B\) being a skew-symmetrizable \((n+m \times n)\) -matrix and with \(\mathbf{x}\) being an \(n\)-tuple of independent elements in the field of rational functions in \(n\) variables.
For the compendium on the cluster algebra and quiver package see Arxiv 1102.4844.
AUTHORS:
- Gregg Musiker: Initial Version
- Christian Stump: Initial Version
- Aram Dermenjian (2015-07-01): Updating ability to not rely solely on clusters
- Jesse Levitt (2015-07-01): Updating ability to not rely solely on clusters
REFERENCES:
[BDP2013] | (1, 2, 3, 4) Thomas Brüstle, Grégoire Dupont, Matthieu Pérotin On Maximal Green Sequences Arxiv 1205.2050 |
[FZ2007] | Sergey Fomin, Andrei Zelevinsky “Cluster Algebras IV: coefficients” Arxiv 0602259 |
See also
For mutation types of cluster seeds, see sage.combinat.cluster_algebra_quiver.quiver_mutation_type.QuiverMutationType()
. Cluster seeds are closely related to sage.combinat.cluster_algebra_quiver.quiver.ClusterQuiver()
.
-
class
sage.combinat.cluster_algebra_quiver.cluster_seed.
ClusterSeed
(data, frozen=None, is_principal=False, user_labels=None, user_labels_prefix='x')¶ Bases:
sage.structure.sage_object.SageObject
The cluster seed associated to an exchange matrix.
INPUT:
data
– can be any of the following:* QuiverMutationType * str - a string representing a QuiverMutationType or a common quiver type (see Examples) * ClusterQuiver * Matrix - a skew-symmetrizable matrix * DiGraph - must be the input data for a quiver * List of edges - must be the edge list of a digraph for a quiver
EXAMPLES:
sage: S = ClusterSeed(['A',5]); S A seed for a cluster algebra of rank 5 of type ['A', 5] sage: S = ClusterSeed(['A',[2,5],1]); S A seed for a cluster algebra of rank 7 of type ['A', [2, 5], 1] sage: T = ClusterSeed( S ); T A seed for a cluster algebra of rank 7 of type ['A', [2, 5], 1] sage: T = ClusterSeed( S._M ); T A seed for a cluster algebra of rank 7 sage: T = ClusterSeed( S.quiver()._digraph ); T A seed for a cluster algebra of rank 7 sage: T = ClusterSeed( S.quiver()._digraph.edges() ); T A seed for a cluster algebra of rank 7 sage: S = ClusterSeed(['B',2]); S A seed for a cluster algebra of rank 2 of type ['B', 2] sage: S = ClusterSeed(['C',2]); S A seed for a cluster algebra of rank 2 of type ['B', 2] sage: S = ClusterSeed(['A', [5,0],1]); S A seed for a cluster algebra of rank 5 of type ['D', 5] sage: S = ClusterSeed(['GR',[3,7]]); S A seed for a cluster algebra of rank 6 of type ['E', 6] sage: S = ClusterSeed(['F', 4, [2,1]]); S A seed for a cluster algebra of rank 6 of type ['F', 4, [1, 2]] sage: S = ClusterSeed(['A',4]); S._use_fpolys True sage: S._use_d_vec True sage: S._use_g_vec True sage: S._use_c_vec True sage: S = ClusterSeed(['A',4]); S.use_fpolys(False); S._use_fpolys False
-
b_matrix
()¶ Returns the \(B\) -matrix of
self
.EXAMPLES:
sage: ClusterSeed(['A',4]).b_matrix() [ 0 1 0 0] [-1 0 -1 0] [ 0 1 0 1] [ 0 0 -1 0] sage: ClusterSeed(['B',4]).b_matrix() [ 0 1 0 0] [-1 0 -1 0] [ 0 1 0 1] [ 0 0 -2 0] sage: ClusterSeed(['D',4]).b_matrix() [ 0 1 0 0] [-1 0 -1 -1] [ 0 1 0 0] [ 0 1 0 0] sage: ClusterSeed(QuiverMutationType([['A',2],['B',2]])).b_matrix() [ 0 1 0 0] [-1 0 0 0] [ 0 0 0 1] [ 0 0 -2 0]
-
b_matrix_class
(depth=+Infinity, up_to_equivalence=True)¶ Returns all \(B\)-matrices in the mutation class of
self
.INPUT:
depth
– (default:infinity) integer or infinity, only seeds with distance at most depth from self are returnedup_to_equivalence
– (default: True) if True, only ‘B’-matrices up to equivalence are considered.
EXAMPLES:
- for examples see
b_matrix_class_iter()
TESTS:
sage: A = ClusterSeed(['A',3]).b_matrix_class() sage: A = ClusterSeed(['A',[2,1],1]).b_matrix_class()
-
b_matrix_class_iter
(depth=+Infinity, up_to_equivalence=True)¶ Returns an iterator through all \(B\)-matrices in the mutation class of
self
.INPUT:
depth
– (default:infinity) integer or infinity, only seeds with distance at most depth from self are returnedup_to_equivalence
– (default: True) if True, only ‘B’-matrices up to equivalence are considered.
EXAMPLES:
A standard finite type example:
sage: S = ClusterSeed(['A',4]) sage: it = S.b_matrix_class_iter() sage: for T in it: print(T) [ 0 0 0 1] [ 0 0 1 1] [ 0 -1 0 0] [-1 -1 0 0] [ 0 0 0 1] [ 0 0 1 0] [ 0 -1 0 1] [-1 0 -1 0] [ 0 0 1 1] [ 0 0 0 -1] [-1 0 0 0] [-1 1 0 0] [ 0 0 0 1] [ 0 0 -1 1] [ 0 1 0 -1] [-1 -1 1 0] [ 0 0 0 1] [ 0 0 -1 0] [ 0 1 0 -1] [-1 0 1 0] [ 0 0 0 -1] [ 0 0 -1 1] [ 0 1 0 -1] [ 1 -1 1 0]
A finite type example with given depth:
sage: it = S.b_matrix_class_iter(depth=1) sage: for T in it: print(T) [ 0 0 0 1] [ 0 0 1 1] [ 0 -1 0 0] [-1 -1 0 0] [ 0 0 0 1] [ 0 0 1 0] [ 0 -1 0 1] [-1 0 -1 0] [ 0 0 1 1] [ 0 0 0 -1] [-1 0 0 0] [-1 1 0 0]
Finite type example not considered up to equivalence:
sage: S = ClusterSeed(['A',3]) sage: it = S.b_matrix_class_iter(up_to_equivalence=False) sage: for T in it: print(T) [ 0 1 0] [-1 0 -1] [ 0 1 0] [ 0 1 0] [-1 0 1] [ 0 -1 0] [ 0 -1 0] [ 1 0 1] [ 0 -1 0] [ 0 -1 0] [ 1 0 -1] [ 0 1 0] [ 0 -1 1] [ 1 0 -1] [-1 1 0] [ 0 1 -1] [-1 0 1] [ 1 -1 0] [ 0 0 1] [ 0 0 -1] [-1 1 0] [ 0 -1 1] [ 1 0 0] [-1 0 0] [ 0 0 -1] [ 0 0 1] [ 1 -1 0] [ 0 1 -1] [-1 0 0] [ 1 0 0] [ 0 1 1] [-1 0 0] [-1 0 0] [ 0 -1 -1] [ 1 0 0] [ 1 0 0] [ 0 0 -1] [ 0 0 -1] [ 1 1 0] [ 0 0 1] [ 0 0 1] [-1 -1 0]
Infinite (but finite mutation) type example:
sage: S = ClusterSeed(['A',[1,2],1]) sage: it = S.b_matrix_class_iter() sage: for T in it: print(T) [ 0 1 1] [-1 0 1] [-1 -1 0] [ 0 -2 1] [ 2 0 -1] [-1 1 0]
Infinite mutation type example:
sage: S = ClusterSeed(['E',10]) sage: it = S.b_matrix_class_iter(depth=3) sage: len ( [T for T in it] ) 266
-
c_matrix
(show_warnings=True)¶ Returns all c-vectors of
self
.Warning: this method assumes the sign-coherence conjecture and that the input seed is sign-coherent (has an exchange matrix with columns of like signs). Otherwise, computational errors might arise.
EXAMPLES:
sage: S = ClusterSeed(['A',3]).principal_extension() sage: S.mutate([2,1,2]) sage: S.c_matrix() [ 1 0 0] [ 0 0 -1] [ 0 -1 0] sage: S = ClusterSeed(['A',4]); sage: S.use_g_vectors(False); S.use_fpolys(False); S.use_c_vectors(False); S.use_d_vectors(False); S.track_mutations(False); sage: S.c_matrix() Traceback (most recent call last): ... ValueError: Unable to calculate c-vectors. Need to use c vectors.
-
c_vector
(k)¶ Returns the
k
-th c-vector ofself
. It is obtained as thek
-th column vector of the bottom part of theB
-matrix ofself
.Warning: this method assumes the sign-coherence conjecture and that the input seed is sign-coherent (has an exchange matrix with columns of like signs). Otherwise, computational errors might arise.
EXAMPLES:
sage: S = ClusterSeed(['A',3]).principal_extension() sage: S.mutate([2,1,2]) sage: [ S.c_vector(k) for k in range(3) ] [(1, 0, 0), (0, 0, -1), (0, -1, 0)] sage: S = ClusterSeed(Matrix([[0,1],[-1,0],[1,0],[-1,1]])); S A seed for a cluster algebra of rank 2 with 2 frozen variables sage: S.c_vector(0) (1, 0) sage: S = ClusterSeed(Matrix([[0,1],[-1,0],[1,0],[-1,1]])); S.use_c_vectors(bot_is_c=True); S A seed for a cluster algebra of rank 2 with 2 frozen variables sage: S.c_vector(0) (1, -1)
-
cluster
()¶ Returns a copy of the cluster of
self
.EXAMPLES:
sage: S = ClusterSeed(['A',3]) sage: S.cluster() [x0, x1, x2] sage: S.mutate(1) sage: S.cluster() [x0, (x0*x2 + 1)/x1, x2] sage: S.mutate(2) sage: S.cluster() [x0, (x0*x2 + 1)/x1, (x0*x2 + x1 + 1)/(x1*x2)] sage: S.mutate([2,1]) sage: S.cluster() [x0, x1, x2]
-
cluster_class
(depth=+Infinity, show_depth=False, up_to_equivalence=True)¶ Returns the cluster class of
self
with respect to certain constraints.INPUT:
depth
– (default: infinity) integer, only seeds with distance at most depth from self are returnedreturn_depth
– (default False) - if True, ignored if depth is set; returns the depth of the mutation class, i.e., the maximal distance from self of an element in the mutation classup_to_equivalence
– (default: True) if True, only clusters up to equivalence are considered.
EXAMPLES:
- for examples see
cluster_class_iter()
TESTS:
sage: A = ClusterSeed(['A',3]).cluster_class()
-
cluster_class_iter
(depth=+Infinity, show_depth=False, up_to_equivalence=True)¶ Returns an iterator through all clusters in the mutation class of
self
.INPUT:
depth
– (default: infinity) integer or infinity, only seeds with distance at most depth from self are returnedshow_depth
– (default False) - if True, ignored if depth is set; returns the depth of the mutation class, i.e., the maximal distance from self of an element in the mutation classup_to_equivalence
– (default: True) if True, only clusters up to equivalence are considered.
EXAMPLES:
A standard finite type example:
sage: S = ClusterSeed(['A',3]) sage: it = S.cluster_class_iter() sage: for T in it: print(T) [x0, x1, x2] [x0, x1, (x1 + 1)/x2] [x0, (x0*x2 + 1)/x1, x2] [(x1 + 1)/x0, x1, x2] [x0, (x0*x2 + x1 + 1)/(x1*x2), (x1 + 1)/x2] [(x1 + 1)/x0, x1, (x1 + 1)/x2] [(x1 + 1)/x0, (x0*x2 + x1 + 1)/(x0*x1), x2] [x0, (x0*x2 + 1)/x1, (x0*x2 + x1 + 1)/(x1*x2)] [(x0*x2 + x1 + 1)/(x0*x1), (x0*x2 + 1)/x1, x2] [(x1 + 1)/x0, (x1^2 + x0*x2 + 2*x1 + 1)/(x0*x1*x2), (x1 + 1)/x2] [(x1 + 1)/x0, (x0*x2 + x1 + 1)/(x0*x1), (x1^2 + x0*x2 + 2*x1 + 1)/(x0*x1*x2)] [(x1^2 + x0*x2 + 2*x1 + 1)/(x0*x1*x2), (x0*x2 + x1 + 1)/(x1*x2), (x1 + 1)/x2] [(x0*x2 + x1 + 1)/(x0*x1), (x0*x2 + 1)/x1, (x0*x2 + x1 + 1)/(x1*x2)] [(x0*x2 + x1 + 1)/(x1*x2), (x0*x2 + x1 + 1)/(x0*x1), (x1^2 + x0*x2 + 2*x1 + 1)/(x0*x1*x2)]
A finite type example with given depth:
sage: it = S.cluster_class_iter(depth=1) sage: for T in it: print(T) [x0, x1, x2] [x0, x1, (x1 + 1)/x2] [x0, (x0*x2 + 1)/x1, x2] [(x1 + 1)/x0, x1, x2]
A finite type example where the depth is returned while computing:
sage: it = S.cluster_class_iter(show_depth=True) sage: for T in it: print(T) [x0, x1, x2] Depth: 0 found: 1 Time: ... s [x0, x1, (x1 + 1)/x2] [x0, (x0*x2 + 1)/x1, x2] [(x1 + 1)/x0, x1, x2] Depth: 1 found: 4 Time: ... s [x0, (x0*x2 + x1 + 1)/(x1*x2), (x1 + 1)/x2] [(x1 + 1)/x0, x1, (x1 + 1)/x2] [(x1 + 1)/x0, (x0*x2 + x1 + 1)/(x0*x1), x2] [x0, (x0*x2 + 1)/x1, (x0*x2 + x1 + 1)/(x1*x2)] [(x0*x2 + x1 + 1)/(x0*x1), (x0*x2 + 1)/x1, x2] Depth: 2 found: 9 Time: ... s [(x1 + 1)/x0, (x1^2 + x0*x2 + 2*x1 + 1)/(x0*x1*x2), (x1 + 1)/x2] [(x1 + 1)/x0, (x0*x2 + x1 + 1)/(x0*x1), (x1^2 + x0*x2 + 2*x1 + 1)/(x0*x1*x2)] [(x1^2 + x0*x2 + 2*x1 + 1)/(x0*x1*x2), (x0*x2 + x1 + 1)/(x1*x2), (x1 + 1)/x2] [(x0*x2 + x1 + 1)/(x0*x1), (x0*x2 + 1)/x1, (x0*x2 + x1 + 1)/(x1*x2)] Depth: 3 found: 13 Time: ... s [(x0*x2 + x1 + 1)/(x1*x2), (x0*x2 + x1 + 1)/(x0*x1), (x1^2 + x0*x2 + 2*x1 + 1)/(x0*x1*x2)] Depth: 4 found: 14 Time: ... s
Finite type examples not considered up to equivalence:
sage: it = S.cluster_class_iter(up_to_equivalence=False) sage: len( [ T for T in it ] ) 84 sage: it = ClusterSeed(['A',2]).cluster_class_iter(up_to_equivalence=False) sage: for T in it: print(T) [x0, x1] [x0, (x0 + 1)/x1] [(x1 + 1)/x0, x1] [(x1 + 1)/x0, (x0 + x1 + 1)/(x0*x1)] [(x0 + x1 + 1)/(x0*x1), (x0 + 1)/x1] [(x0 + x1 + 1)/(x0*x1), (x1 + 1)/x0] [(x0 + 1)/x1, (x0 + x1 + 1)/(x0*x1)] [x1, (x1 + 1)/x0] [(x0 + 1)/x1, x0] [x1, x0]
Infinite type examples:
sage: S = ClusterSeed(['A',[1,1],1]) sage: it = S.cluster_class_iter() sage: next(it) [x0, x1] sage: next(it) [x0, (x0^2 + 1)/x1] sage: next(it) [(x1^2 + 1)/x0, x1] sage: next(it) [(x0^4 + 2*x0^2 + x1^2 + 1)/(x0*x1^2), (x0^2 + 1)/x1] sage: next(it) [(x1^2 + 1)/x0, (x1^4 + x0^2 + 2*x1^2 + 1)/(x0^2*x1)] sage: it = S.cluster_class_iter(depth=3) sage: for T in it: print(T) [x0, x1] [x0, (x0^2 + 1)/x1] [(x1^2 + 1)/x0, x1] [(x0^4 + 2*x0^2 + x1^2 + 1)/(x0*x1^2), (x0^2 + 1)/x1] [(x1^2 + 1)/x0, (x1^4 + x0^2 + 2*x1^2 + 1)/(x0^2*x1)] [(x0^4 + 2*x0^2 + x1^2 + 1)/(x0*x1^2), (x0^6 + 3*x0^4 + 2*x0^2*x1^2 + x1^4 + 3*x0^2 + 2*x1^2 + 1)/(x0^2*x1^3)] [(x1^6 + x0^4 + 2*x0^2*x1^2 + 3*x1^4 + 2*x0^2 + 3*x1^2 + 1)/(x0^3*x1^2), (x1^4 + x0^2 + 2*x1^2 + 1)/(x0^2*x1)]
-
cluster_index
(cluster_str)¶ Returns the index of a cluster if use_fpolys is on
INPUT:
cluster_str
– The string to look for in the cluster
OUTPUT:
Returns an integer or None if the string is not a cluster variable
EXAMPLES:
sage: S = ClusterSeed(['A',4],user_labels=['x','y','z','w']); S.mutate('x') sage: S.cluster_index('x') sage: S.cluster_index('(y+1)/x') 0
-
cluster_variable
(k)¶ Generates a cluster variable using F-polynomials
EXAMPLES:
sage: S = ClusterSeed(['A',3]) sage: S.mutate([0,1]) sage: S.cluster_variable(0) (x1 + 1)/x0 sage: S.cluster_variable(1) (x0*x2 + x1 + 1)/(x0*x1)
-
coefficient
(k)¶ Returns the coefficient of
self
at indexk
.EXAMPLES:
sage: S = ClusterSeed(['A',3]).principal_extension() sage: S.mutate([2,1,2]) sage: [ S.coefficient(k) for k in range(3) ] [y0, 1/y2, 1/y1]
-
coefficients
()¶ Returns all coefficients of
self
.EXAMPLES:
sage: S = ClusterSeed(['A',3]).principal_extension() sage: S.mutate([2,1,2]) sage: S.coefficients() [y0, 1/y2, 1/y1]
-
d_matrix
(show_warnings=True)¶ Returns the matrix of d-vectors of
self
.- EXAMPLES::
- sage: S = ClusterSeed([‘A’,4]); S.d_matrix() [-1 0 0 0] [ 0 -1 0 0] [ 0 0 -1 0] [ 0 0 0 -1] sage: S.mutate([1,2,1,0,1,3]); S.d_matrix() [1 1 0 1] [1 1 1 1] [1 0 1 1] [0 0 0 1]
-
d_vector
(k)¶ Returns the
k
-th d-vector ofself
. This is the exponent vector of the denominator of thek
-th cluster variable.EXAMPLES:
sage: S = ClusterSeed(['A',3]) sage: S.mutate([2,1,2]) sage: [ S.d_vector(k) for k in range(3) ] [(-1, 0, 0), (0, 1, 1), (0, 1, 0)]
-
exchangeable_part
()¶ Returns the restriction to the principal part (i.e. the exchangeable variables) of
self
.EXAMPLES:
sage: S = ClusterSeed(['A',4]) sage: T = ClusterSeed( S.quiver().digraph().edges(), frozen=1 ) sage: T.quiver().digraph().edges() [(0, 1, (1, -1)), (2, 1, (1, -1)), (2, 3, (1, -1))] sage: T.exchangeable_part().quiver().digraph().edges() [(0, 1, (1, -1)), (2, 1, (1, -1))]
-
f_polynomial
(k)¶ Returns the
k
-th F-polynomial ofself
. It is obtained from thek
-th cluster variable by setting all \(x_i\) to \(1\).Warning: this method assumes the sign-coherence conjecture and that the input seed is sign-coherent (has an exchange matrix with columns of like signs). Otherwise, computational errors might arise.
EXAMPLES:
sage: S = ClusterSeed(['A',3]).principal_extension() sage: S.mutate([2,1,2]) sage: [S.f_polynomial(k) for k in range(3)] [1, y1*y2 + y2 + 1, y1 + 1] sage: S = ClusterSeed(Matrix([[0,1],[-1,0],[1,0],[-1,1]])); S.use_c_vectors(bot_is_c=True); S A seed for a cluster algebra of rank 2 with 2 frozen variables sage: T = ClusterSeed(Matrix([[0,1],[-1,0]])).principal_extension(); T A seed for a cluster algebra of rank 2 with principal coefficients sage: S.mutate(0) sage: T.mutate(0) sage: S.f_polynomials() [y0 + y1, 1] sage: T.f_polynomials() [y0 + 1, 1]
-
f_polynomials
()¶ Returns all F-polynomials of
self
. These are obtained from the cluster variables by setting all \(x_i\)‘s to \(1\).Warning: this method assumes the sign-coherence conjecture and that the input seed is sign-coherent (has an exchange matrix with columns of like signs). Otherwise, computational errors might arise.
EXAMPLES:
sage: S = ClusterSeed(['A',3]).principal_extension() sage: S.mutate([2,1,2]) sage: S.f_polynomials() [1, y1*y2 + y2 + 1, y1 + 1]
-
find_upper_bound
(verbose=False)¶ Return the upper bound of the given cluster algebra as a quotient_ring.
The upper bound is the intersection of the Laurent polynomial rings of the initial cluster and its neighboring clusters. As such, it always contains both the cluster algebra and the upper cluster algebra. This function uses the algorithm from [MaMu].
When the initial seed is totally coprime (for example, when the unfrozen part of the exchange matrix has full rank), the upper bound is equal to the upper cluster algebra by [BeFoZe].
Warning
The computation time grows rapidly with the size of the seed and the number of steps. For most seeds larger than four vertices, the algorithm may take an infeasible amount of time. Additionally, it will run forever without terminating whenever the upper bound is infinitely-generated (such as the example in [Sp]).
INPUT:
verbose
– (default:False
) ifTrue
, prints output during the computation.
EXAMPLES:
finite type:
sage: S = ClusterSeed(['A',3]) sage: S.find_upper_bound() Quotient of Multivariate Polynomial Ring in x0, x1, x2, x0p, x1p, x2p, z0 over Rational Field by the ideal (x0*x0p - x1 - 1, x1*x1p - x0*x2 - 1, x2*x2p - x1 - 1, x0*z0 - x2p, x1*z0 + z0 - x0p*x2p, x2*z0 - x0p, x1p*z0 + z0 - x0p*x1p*x2p + x1 + 1)
Markov:
sage: B = matrix([[0,2,-2],[-2,0,2],[2,-2,0]]) sage: S = ClusterSeed(B) sage: S.find_upper_bound() Quotient of Multivariate Polynomial Ring in x0, x1, x2, x0p, x1p, x2p, z0 over Rational Field by the ideal (x0*x0p - x2^2 - x1^2, x1*x1p - x2^2 - x0^2, x2*x2p - x1^2 - x0^2, x0p*x1p*x2p - x0*x1*x2p - x0*x2*x1p - x1*x2*x0p - 2*x0*x1*x2, x0^3*z0 - x1p*x2p + x1*x2, x0*x1*z0 - x2p - x2, x1^3*z0 - x0p*x2p + x0*x2, x0*x2*z0 - x1p - x1, x1*x2*z0 - x0p - x0, x2^3*z0 - x0p*x1p + x0*x1)
REFERENCES:
[BeFoZe] Berenstein-Fomin-Zelevinsky, Cluster algebras III: Upper bounds and double Bruhat cells, Duke Math. J., 126(1):1-52, 2005. [MaMu] Matherne-Muller, Computing upper cluster algebras, Int. Math. Res. Not. IMRN, 2015, 3121-3149. [Sp] Speyer, An infinitely generated upper cluster algebra, Arxiv 1305.6867.
-
first_green_vertex
()¶ Return the first green vertex of
self
.A vertex is defined to be green if its c-vector has all non-positive entries. More information on green vertices can be found at [BDP2013]
EXAMPLES:
sage: ClusterSeed(['A',3]).principal_extension().first_green_vertex() 0 sage: ClusterSeed(['A',[3,3],1]).principal_extension().first_green_vertex() 0
-
first_red_vertex
()¶ Return the first red vertex of
self
.A vertex is defined to be red if its c-vector has all non-negative entries. More information on red vertices can be found at [BDP2013].
EXAMPLES:
sage: ClusterSeed(['A',3]).principal_extension().first_red_vertex() sage: ClusterSeed(['A',[3,3],1]).principal_extension().first_red_vertex() sage: Q = ClusterSeed(['A',[3,3],1]).principal_extension(); sage: Q.mutate(1); sage: Q.first_red_vertex() 1
-
first_urban_renewal
()¶ Return the first urban renewal vertex.
An urban renewal vertex is one in which there are two arrows pointing toward the vertex and two arrows pointing away.
EXAMPLES:
sage: G = ClusterSeed(['GR',[4,9]]); G.first_urban_renewal() 5
-
g_matrix
(show_warnings=True)¶ Returns the matrix of all g-vectors of
self
. These are the degree vectors of the cluster variables after setting all \(y_i\)‘s to \(0\).Warning: this method assumes the sign-coherence conjecture and that the input seed is sign-coherent (has an exchange matrix with columns of like signs). Otherwise, computational errors might arise.
EXAMPLES:
sage: S = ClusterSeed(['A',3]).principal_extension() sage: S.mutate([2,1,2]) sage: S.g_matrix() [ 1 0 0] [ 0 0 -1] [ 0 -1 0] sage: S = ClusterSeed(['A',3]) sage: S.mutate([0,1]) sage: S.g_matrix() [-1 -1 0] [ 1 0 0] [ 0 0 1] sage: S = ClusterSeed(['A',4]); S.use_g_vectors(False); S.use_fpolys(False); S.g_matrix() [1 0 0 0] [0 1 0 0] [0 0 1 0] [0 0 0 1] sage: S = ClusterSeed(['A',4]) sage: S.use_g_vectors(False); S.use_c_vectors(False); S.use_fpolys(False); S.track_mutations(False); S.g_matrix() Traceback (most recent call last): ... ValueError: Unable to calculate g-vectors. Need to use g vectors.
-
g_vector
(k)¶ Returns the
k
-th g-vector ofself
. This is the degree vector of thek
-th cluster variable after setting all \(y_i\)‘s to \(0\).Warning: this method assumes the sign-coherence conjecture and that the input seed is sign-coherent (has an exchange matrix with columns of like signs). Otherwise, computational errors might arise.
EXAMPLES:
sage: S = ClusterSeed(['A',3]).principal_extension() sage: S.mutate([2,1,2]) sage: [ S.g_vector(k) for k in range(3) ] [(1, 0, 0), (0, 0, -1), (0, -1, 0)]
-
greedy
(a1, a2, algorithm='by_recursion')¶ Returns the greedy element \(x[a_1,a_2]\) assuming that self is rank two.
The third input can be ‘by_recursion’, ‘by_combinatorics’, or ‘just_numbers’ to specify if the user wants the element computed by the recurrence, combinatorial formula, or wants to set \(x_1\) and \(x_2\) to be one.
See [LeeLiZe] for more details.
EXAMPLES:
sage: S = ClusterSeed(['R2', [3, 3]]) sage: S.greedy(4, 4) (x0^12 + x1^12 + 4*x0^9 + 4*x1^9 + 6*x0^6 + 4*x0^3*x1^3 + 6*x1^6 + 4*x0^3 + 4*x1^3 + 1)/(x0^4*x1^4) sage: S.greedy(4, 4, 'by_combinatorics') (x0^12 + x1^12 + 4*x0^9 + 4*x1^9 + 6*x0^6 + 4*x0^3*x1^3 + 6*x1^6 + 4*x0^3 + 4*x1^3 + 1)/(x0^4*x1^4) sage: S.greedy(4, 4, 'just_numbers') 35 sage: S = ClusterSeed(['R2', [2, 2]]) sage: S.greedy(1, 2) (x0^4 + 2*x0^2 + x1^2 + 1)/(x0*x1^2) sage: S.greedy(1, 2, 'by_combinatorics') (x0^4 + 2*x0^2 + x1^2 + 1)/(x0*x1^2)
REFERENCES:
[LeeLiZe] (1, 2) Lee-Li-Zelevinsky, Greedy elements in rank 2 cluster algebras, Arxiv 1208.2391
-
green_vertices
()¶ Return the list of green vertices of
self
.A vertex is defined to be green if its c-vector has all non-positive entries. More information on green vertices can be found at [BDP2013]
OUTPUT:
The green vertices as a list of integers.
EXAMPLES:
sage: ClusterSeed(['A',3]).principal_extension().green_vertices() [0, 1, 2] sage: ClusterSeed(['A',[3,3],1]).principal_extension().green_vertices() [0, 1, 2, 3, 4, 5]
-
ground_field
()¶ Returns the ground field of the cluster of
self
.EXAMPLES:
sage: S = ClusterSeed(['A',3]) sage: S.ground_field() Multivariate Polynomial Ring in x0, x1, x2, y0, y1, y2 over Rational Field
-
highest_degree_denominator
(filter=None)¶ Return the vertex of the cluster polynomial with highest degree in the denominator.
INPUT:
filter
- Filter should be a list or iterable
OUTPUT:
Returns an integer.
EXAMPLES:
sage: B = matrix([[0,-1,0,-1,1,1],[1,0,1,0,-1,-1],[0,-1,0,-1,1,1],[1,0,1,0,-1,-1],[-1,1,-1,1,0,0],[-1,1,-1,1,0,0]]) sage: C = ClusterSeed(B).principal_extension(); C.mutate([0,1,2,4,3,2,5,4,3]) sage: C.highest_degree_denominator() 5
-
interact
(fig_size=1, circular=True)¶ Only in notebook mode. Starts an interactive window for cluster seed mutations.
INPUT:
fig_size
– (default: 1) factor by which the size of the plot is multiplied.circular
– (default: True) if True, the circular plot is chosen, otherwise >>spring<< is used.
TESTS:
sage: S = ClusterSeed(['A',4]) sage: S.interact() # long time 'The interactive mode only runs in the Sage notebook.'
-
is_acyclic
()¶ Returns True iff self is acyclic (i.e., if the underlying quiver is acyclic).
EXAMPLES:
sage: ClusterSeed(['A',4]).is_acyclic() True sage: ClusterSeed(['A',[2,1],1]).is_acyclic() True sage: ClusterSeed([[0,1],[1,2],[2,0]]).is_acyclic() False
-
is_bipartite
(return_bipartition=False)¶ Returns True iff self is bipartite (i.e., if the underlying quiver is bipartite).
INPUT:
- return_bipartition – (default:False) if True, the bipartition is returned in the case of
self
being bipartite.
EXAMPLES:
sage: ClusterSeed(['A',[3,3],1]).is_bipartite() True sage: ClusterSeed(['A',[4,3],1]).is_bipartite() False
- return_bipartition – (default:False) if True, the bipartition is returned in the case of
-
is_finite
()¶ Returns True if
self
is of finite type.EXAMPLES:
sage: S = ClusterSeed(['A',3]) sage: S.is_finite() True sage: S = ClusterSeed(['A',[2,2],1]) sage: S.is_finite() False
-
is_mutation_finite
(nr_of_checks=None, return_path=False)¶ Returns True if
self
is of finite mutation type.INPUT:
nr_of_checks
– (default: None) number of mutations applied. Standard is 500*(number of vertices of self).return_path
– (default: False) if True, in case of self not being mutation finite, a path from self to a quiver with an edge label (a,-b) and a*b > 4 is returned.
ALGORITHM:
- A cluster seed is mutation infinite if and only if every \(b_{ij}*b_{ji} > -4\). Thus, we apply random mutations in random directions
WARNING:
- Uses a non-deterministic method by random mutations in various directions.
- In theory, it can return a wrong True.
EXAMPLES:
sage: S = ClusterSeed(['A',10]) sage: S._mutation_type = None sage: S.is_mutation_finite() True sage: S = ClusterSeed([(0,1),(1,2),(2,3),(3,4),(4,5),(5,6),(6,7),(7,8),(2,9)]) sage: S.is_mutation_finite() False
-
m
()¶ Returns the number of frozen variables of
self
.EXAMPLES:
sage: S = ClusterSeed(['A',3]) sage: S.n() 3 sage: S.m() 0 sage: S = S.principal_extension() sage: S.m() 3
-
most_decreased_denominator_after_mutation
()¶ Return the vertex that will produce the most decrease in denominator degrees after mutation
EXAMPLES:
sage: S = ClusterSeed(['A',5]) sage: S.mutate([0,2,3,1,2,3,1,2,0,2,3]) sage: S.most_decreased_denominator_after_mutation() 2
-
most_decreased_edge_after_mutation
()¶ Return the vertex that will produce the least degrees after mutation
EXAMPLES:
sage: S = ClusterSeed(['A',5]) sage: S.mutate([0,2,3,1,2,3,1,2,0,2,3]) sage: S.most_decreased_edge_after_mutation() 2
-
mutate
(sequence, inplace=True)¶ Mutates
self
at a vertex or a sequence of vertices.INPUT:
sequence
– a vertex ofself
, an iterator of vertices ofself
, a function which takes in the ClusterSeed and returns a vertex or an iterator of vertices, or a string representing a type of vertices to mutate.inplace
– (default: True) if False, the result is returned, otherwiseself
is modified.
Possible values for vertex types in
sequence
are:"first_source"
: mutates at first found source vertex,"sources"
: mutates at all sources,"first_sink"
: mutates at first sink,"sinks"
: mutates at all sink vertices,"green"
: mutates at the first green vertex,"red"
: mutates at the first red vertex,"urban_renewal"
or"urban"
: mutates at first urban renewal vertex,"all_urban_renewals"
or"all_urban"
: mutates at all urban renewal vertices.
EXAMPLES:
sage: S = ClusterSeed(['A',4]); S.b_matrix() [ 0 1 0 0] [-1 0 -1 0] [ 0 1 0 1] [ 0 0 -1 0] sage: S.mutate(0); S.b_matrix() [ 0 -1 0 0] [ 1 0 -1 0] [ 0 1 0 1] [ 0 0 -1 0] sage: T = S.mutate(0, inplace=False); T A seed for a cluster algebra of rank 4 of type ['A', 4] sage: S.mutate(0) sage: S == T True sage: S.mutate([0,1,0]) sage: S.b_matrix() [ 0 -1 1 0] [ 1 0 0 0] [-1 0 0 1] [ 0 0 -1 0] sage: S = ClusterSeed(QuiverMutationType([['A',1],['A',3]])) sage: S.b_matrix() [ 0 0 0 0] [ 0 0 1 0] [ 0 -1 0 -1] [ 0 0 1 0] sage: T = S.mutate(0,inplace=False) sage: S == T False sage: Q = ClusterSeed(['A',3]);Q.b_matrix() [ 0 1 0] [-1 0 -1] [ 0 1 0] sage: Q.mutate('first_sink');Q.b_matrix() [ 0 -1 0] [ 1 0 1] [ 0 -1 0] sage: def last_vertex(self): return self._n - 1 sage: Q.mutate(last_vertex); Q.b_matrix() [ 0 -1 0] [ 1 0 -1] [ 0 1 0] sage: S = ClusterSeed(['A',4], user_labels=['a','b','c','d']); sage: S.mutate('a'); S.mutate('(b+1)/a') sage: S.cluster() [a, b, c, d] sage: S = ClusterSeed(['A',4], user_labels=['a','b','c']); Traceback (most recent call last): ... ValueError: The number of user-defined labels is not the number of exchangeable and frozen variables. sage: S = ClusterSeed(['A',4],user_labels=['x','y','w','z']) sage: S.mutate('x') sage: S.cluster() [(y + 1)/x, y, w, z] sage: S.mutate('(y+1)/x') sage: S.cluster() [x, y, w, z] sage: S.mutate('y') sage: S.cluster() [x, (x*w + 1)/y, w, z] sage: S.mutate('(x*w+1)/y') sage: S.cluster() [x, y, w, z] sage: S = ClusterSeed(['A',4], user_labels=[[1,2],[2,3],[4,5],[5,6]]); sage: S.cluster() [x_1_2, x_2_3, x_4_5, x_5_6] sage: S.mutate('[1,2]'); sage: S.cluster() [(x_2_3 + 1)/x_1_2, x_2_3, x_4_5, x_5_6] sage: S = ClusterSeed(['A',4], user_labels=[[1,2],[2,3],[4,5],[5,6]],user_labels_prefix='P'); sage: S.cluster() [P_1_2, P_2_3, P_4_5, P_5_6] sage: S.mutate('[1,2]') sage: S.cluster() [(P_2_3 + 1)/P_1_2, P_2_3, P_4_5, P_5_6] sage: S.mutate('P_4_5') sage: S.cluster() [(P_2_3 + 1)/P_1_2, P_2_3, (P_2_3*P_5_6 + 1)/P_4_5, P_5_6] sage: S = ClusterSeed(['A',4]) sage: S.mutate([0,1,0,1,0,2,1]) sage: T = ClusterSeed(S) sage: S.use_fpolys(False) sage: S.use_g_vectors(False) sage: S.use_c_vectors(False) sage: S._C sage: S._G sage: S._F sage: S.g_matrix() [ 0 -1 0 0] [ 1 1 1 0] [ 0 0 -1 0] [ 0 0 1 1] sage: S.c_matrix() [ 1 -1 0 0] [ 1 0 0 0] [ 1 0 -1 1] [ 0 0 0 1] sage: S.f_polynomials() == T.f_polynomials() True sage: S.cluster() == T.cluster() True sage: S._mut_path [0, 1, 0, 1, 0, 2, 1]
-
mutation_analysis
(options=['all'], filter=None)¶ Runs an analysis of all potential mutation options. Note that this might take a long time on large seeds.
Notes: Edges are only returned if we have a non-valued quiver. Green and red vertices are only returned if the cluster is principal.
INPUT:
options
– (default: [‘all’]) a list of mutation options.filter
– (default: None) A vertex or interval of vertices to limit our search to
Possible options are:
"all"
- All options below"edges"
- Number of edges (works with skew-symmetric quivers)"edge_diff"
- Edges added/deleted (works with skew-symmetric quivers)"green_vertices"
- List of green vertices (works with principals)"green_vertices_diff"
- Green vertices added/removed (works with principals)"red_vertices"
- List of red vertices (works with principals)"red_vertices_diff"
- Red vertices added/removed (works with principals)"urban_renewals"
- List of urban renewal vertices"urban_renewals_diff"
- Urban renewal vertices added/removed"sources"
- List of source vertices"sources_diff"
- Source vertices added/removed"sinks"
- List of sink vertices"sinks_diff"
- Sink vertices added/removed"denominators"
- List of all denominators of the cluster variables
OUTPUT:
Outputs a dictionary indexed by the vertex numbers. Each vertex will itself also be a dictionary with each desired option included as a key in the dictionary. As an example you would get something similar to: {0: {‘edges’: 1}, 1: {‘edges’: 2}}. This represents that if you were to do a mutation at the current seed then mutating at vertex 0 would result in a quiver with 1 edge and mutating at vertex 0 would result in a quiver with 2 edges.
EXAMPLES:
sage: B = [[0, 4, 0, -1],[-4,0, 3, 0],[0, -3, 0, 1],[1, 0, -1, 0]] sage: S = ClusterSeed(matrix(B)); S.mutate([2,3,1,2,1,3,0,2]) sage: S.mutation_analysis() {0: {'d_matrix': [ 0 0 1 0] [ 0 -1 0 0] [ 0 0 0 -1] [-1 0 0 0], 'denominators': [1, 1, x0, 1], 'edge_diff': 6, 'edges': 13, 'green_vertices': [0, 1, 3], 'green_vertices_diff': {'added': [0], 'removed': []}, 'red_vertices': [2], 'red_vertices_diff': {'added': [], 'removed': [0]}, 'sinks': [], 'sinks_diff': {'added': [], 'removed': [2]}, 'sources': [], 'sources_diff': {'added': [], 'removed': []}, 'urban_renewals': [], 'urban_renewals_diff': {'added': [], 'removed': []}}, 1: {'d_matrix': [ 1 4 1 0] [ 0 1 0 0] [ 0 0 0 -1] [ 1 4 0 0], 'denominators': [x0*x3, x0^4*x1*x3^4, x0, 1], 'edge_diff': 2, 'edges': 9, 'green_vertices': [0, 3], 'green_vertices_diff': {'added': [0], 'removed': [1]}, 'red_vertices': [1, 2], 'red_vertices_diff': {'added': [1], 'removed': [0]}, 'sinks': [2], 'sinks_diff': {'added': [], 'removed': []}, 'sources': [], 'sources_diff': {'added': [], 'removed': []}, 'urban_renewals': [], 'urban_renewals_diff': {'added': [], 'removed': []}}, 2: {'d_matrix': [ 1 0 0 0] [ 0 -1 0 0] [ 0 0 0 -1] [ 1 0 1 0], 'denominators': [x0*x3, 1, x3, 1], 'edge_diff': 0, 'edges': 7, 'green_vertices': [1, 2, 3], 'green_vertices_diff': {'added': [2], 'removed': []}, 'red_vertices': [0], 'red_vertices_diff': {'added': [], 'removed': [2]}, 'sinks': [], 'sinks_diff': {'added': [], 'removed': [2]}, 'sources': [2], 'sources_diff': {'added': [2], 'removed': []}, 'urban_renewals': [], 'urban_renewals_diff': {'added': [], 'removed': []}}, 3: {'d_matrix': [ 1 0 1 1] [ 0 -1 0 0] [ 0 0 0 1] [ 1 0 0 1], 'denominators': [x0*x3, 1, x0, x0*x2*x3], 'edge_diff': -1, 'edges': 6, 'green_vertices': [1], 'green_vertices_diff': {'added': [], 'removed': [3]}, 'red_vertices': [0, 2, 3], 'red_vertices_diff': {'added': [3], 'removed': []}, 'sinks': [2], 'sinks_diff': {'added': [], 'removed': []}, 'sources': [1], 'sources_diff': {'added': [1], 'removed': []}, 'urban_renewals': [], 'urban_renewals_diff': {'added': [], 'removed': []}}} sage: S = ClusterSeed(['A',3]).principal_extension() sage: S.mutation_analysis() {0: {'d_matrix': [ 1 0 0] [ 0 -1 0] [ 0 0 -1], 'denominators': [x0, 1, 1], 'green_vertices': [1, 2], 'green_vertices_diff': {'added': [], 'removed': [0]}, 'red_vertices': [0], 'red_vertices_diff': {'added': [0], 'removed': []}, 'sinks': [], 'sinks_diff': {'added': [], 'removed': [1]}, 'sources': [4, 5], 'sources_diff': {'added': [], 'removed': [3]}, 'urban_renewals': [], 'urban_renewals_diff': {'added': [], 'removed': []}}, 1: {'d_matrix': [-1 0 0] [ 0 1 0] [ 0 0 -1], 'denominators': [1, x1, 1], 'green_vertices': [0, 2], 'green_vertices_diff': {'added': [], 'removed': [1]}, 'red_vertices': [1], 'red_vertices_diff': {'added': [1], 'removed': []}, 'sinks': [0, 2, 4], 'sinks_diff': {'added': [0, 2, 4], 'removed': [1]}, 'sources': [1, 3, 5], 'sources_diff': {'added': [1], 'removed': [4]}, 'urban_renewals': [], 'urban_renewals_diff': {'added': [], 'removed': []}}, 2: {'d_matrix': [-1 0 0] [ 0 -1 0] [ 0 0 1], 'denominators': [1, 1, x2], 'green_vertices': [0, 1], 'green_vertices_diff': {'added': [], 'removed': [2]}, 'red_vertices': [2], 'red_vertices_diff': {'added': [2], 'removed': []}, 'sinks': [], 'sinks_diff': {'added': [], 'removed': [1]}, 'sources': [3, 4], 'sources_diff': {'added': [], 'removed': [5]}, 'urban_renewals': [], 'urban_renewals_diff': {'added': [], 'removed': []}}}
-
mutation_class
(depth=+Infinity, show_depth=False, return_paths=False, up_to_equivalence=True, only_sink_source=False)¶ Returns the mutation class of
self
with respect to certain constraints.INPUT:
depth
– (default: infinity) integer, only seeds with distance at most depth from self are returned.show_depth
– (default: False) if True, the actual depth of the mutation is shown.return_paths
– (default: False) if True, a shortest path of mutation sequences from self to the given quiver is returned as well.up_to_equivalence
– (default: True) if True, only seeds up to equivalence are considered.sink_source
– (default: False) if True, only mutations at sinks and sources are applied.
EXAMPLES:
- for examples see
mutation_class_iter()
TESTS:
sage: A = ClusterSeed(['A',3]).mutation_class()
-
mutation_class_iter
(depth=+Infinity, show_depth=False, return_paths=False, up_to_equivalence=True, only_sink_source=False)¶ Returns an iterator for the mutation class of
self
with respect to certain constrains.INPUT:
depth
– (default: infinity) integer or infinity, only seeds with distance at mostdepth
fromself
are returned.show_depth
– (default: False) if True, the current depth of the mutation is shown while computing.return_paths
– (default: False) if True, a shortest path of mutations fromself
to the given quiver is returned as well.up_to_equivalence
– (default: True) if True, only one seed up to simultaneous permutation of rows and columns of the exchange matrix is recorded.sink_source
– (default: False) if True, only mutations at sinks and sources are applied.
EXAMPLES:
A standard finite type example:
sage: S = ClusterSeed(['A',3]) sage: it = S.mutation_class_iter() sage: for T in it: print(T) A seed for a cluster algebra of rank 3 of type ['A', 3] A seed for a cluster algebra of rank 3 of type ['A', 3] A seed for a cluster algebra of rank 3 of type ['A', 3] A seed for a cluster algebra of rank 3 of type ['A', 3] A seed for a cluster algebra of rank 3 of type ['A', 3] A seed for a cluster algebra of rank 3 of type ['A', 3] A seed for a cluster algebra of rank 3 of type ['A', 3] A seed for a cluster algebra of rank 3 of type ['A', 3] A seed for a cluster algebra of rank 3 of type ['A', 3] A seed for a cluster algebra of rank 3 of type ['A', 3] A seed for a cluster algebra of rank 3 of type ['A', 3] A seed for a cluster algebra of rank 3 of type ['A', 3] A seed for a cluster algebra of rank 3 of type ['A', 3] A seed for a cluster algebra of rank 3 of type ['A', 3]
A finite type example with given depth:
sage: it = S.mutation_class_iter(depth=1) sage: for T in it: print(T) A seed for a cluster algebra of rank 3 of type ['A', 3] A seed for a cluster algebra of rank 3 of type ['A', 3] A seed for a cluster algebra of rank 3 of type ['A', 3] A seed for a cluster algebra of rank 3 of type ['A', 3]
A finite type example where the depth is shown while computing:
sage: it = S.mutation_class_iter(show_depth=True) sage: for T in it: pass Depth: 0 found: 1 Time: ... s Depth: 1 found: 4 Time: ... s Depth: 2 found: 9 Time: ... s Depth: 3 found: 13 Time: ... s Depth: 4 found: 14 Time: ... s
A finite type example with shortest paths returned:
sage: it = S.mutation_class_iter(return_paths=True) sage: for T in it: print(T) (A seed for a cluster algebra of rank 3 of type ['A', 3], []) (A seed for a cluster algebra of rank 3 of type ['A', 3], [2]) (A seed for a cluster algebra of rank 3 of type ['A', 3], [1]) (A seed for a cluster algebra of rank 3 of type ['A', 3], [0]) (A seed for a cluster algebra of rank 3 of type ['A', 3], [2, 1]) (A seed for a cluster algebra of rank 3 of type ['A', 3], [0, 2]) (A seed for a cluster algebra of rank 3 of type ['A', 3], [0, 1]) (A seed for a cluster algebra of rank 3 of type ['A', 3], [1, 2]) (A seed for a cluster algebra of rank 3 of type ['A', 3], [1, 0]) (A seed for a cluster algebra of rank 3 of type ['A', 3], [0, 2, 1]) (A seed for a cluster algebra of rank 3 of type ['A', 3], [0, 1, 2]) (A seed for a cluster algebra of rank 3 of type ['A', 3], [2, 1, 0]) (A seed for a cluster algebra of rank 3 of type ['A', 3], [1, 0, 2]) (A seed for a cluster algebra of rank 3 of type ['A', 3], [0, 1, 2, 0])
Finite type examples not considered up to equivalence:
sage: it = S.mutation_class_iter(up_to_equivalence=False) sage: len( [ T for T in it ] ) 84 sage: it = ClusterSeed(['A',2]).mutation_class_iter(return_paths=True,up_to_equivalence=False) sage: for T in it: print(T) (A seed for a cluster algebra of rank 2 of type ['A', 2], []) (A seed for a cluster algebra of rank 2 of type ['A', 2], [1]) (A seed for a cluster algebra of rank 2 of type ['A', 2], [0]) (A seed for a cluster algebra of rank 2 of type ['A', 2], [0, 1]) (A seed for a cluster algebra of rank 2 of type ['A', 2], [1, 0]) (A seed for a cluster algebra of rank 2 of type ['A', 2], [1, 0, 1]) (A seed for a cluster algebra of rank 2 of type ['A', 2], [0, 1, 0]) (A seed for a cluster algebra of rank 2 of type ['A', 2], [1, 0, 1, 0]) (A seed for a cluster algebra of rank 2 of type ['A', 2], [0, 1, 0, 1]) (A seed for a cluster algebra of rank 2 of type ['A', 2], [1, 0, 1, 0, 1])
Check that trac ticket #14638 is fixed:
sage: S = ClusterSeed(['E',6]) sage: MC = S.mutation_class(depth=7); len(MC) 534
Infinite type examples:
sage: S = ClusterSeed(['A',[1,1],1]) sage: it = S.mutation_class_iter() sage: next(it) A seed for a cluster algebra of rank 2 of type ['A', [1, 1], 1] sage: next(it) A seed for a cluster algebra of rank 2 of type ['A', [1, 1], 1] sage: next(it) A seed for a cluster algebra of rank 2 of type ['A', [1, 1], 1] sage: next(it) A seed for a cluster algebra of rank 2 of type ['A', [1, 1], 1] sage: it = S.mutation_class_iter(depth=3, return_paths=True) sage: for T in it: print(T) (A seed for a cluster algebra of rank 2 of type ['A', [1, 1], 1], []) (A seed for a cluster algebra of rank 2 of type ['A', [1, 1], 1], [1]) (A seed for a cluster algebra of rank 2 of type ['A', [1, 1], 1], [0]) (A seed for a cluster algebra of rank 2 of type ['A', [1, 1], 1], [1, 0]) (A seed for a cluster algebra of rank 2 of type ['A', [1, 1], 1], [0, 1]) (A seed for a cluster algebra of rank 2 of type ['A', [1, 1], 1], [1, 0, 1]) (A seed for a cluster algebra of rank 2 of type ['A', [1, 1], 1], [0, 1, 0])
-
mutation_sequence
(sequence, show_sequence=False, fig_size=1.2, return_output='seed')¶ Returns the seeds obtained by mutating
self
at all vertices insequence
.INPUT:
sequence
– an iterable of vertices of self.show_sequence
– (default: False) if True, a png containing the associated quivers is shown.fig_size
– (default: 1.2) factor by which the size of the plot is multiplied.return_output
– (default: ‘seed’) determines what output is to be returned:* if 'seed', outputs all the cluster seeds obtained by the ``sequence`` of mutations. * if 'matrix', outputs a list of exchange matrices. * if 'var', outputs a list of new cluster variables obtained at each step.
EXAMPLES:
sage: S = ClusterSeed(['A',2]) sage: for T in S.mutation_sequence([0,1,0]): ....: print(T.b_matrix()) [ 0 -1] [ 1 0] [ 0 1] [-1 0] [ 0 -1] [ 1 0] sage: S=ClusterSeed(['A',2]) sage: S.mutation_sequence([0,1,0,1],return_output='var') [(x1 + 1)/x0, (x0 + x1 + 1)/(x0*x1), (x0 + 1)/x1, x0]
-
mutation_type
()¶ Returns the mutation_type of each connected component of
self
, if it can be determined. Otherwise, the mutation type of this component is set to be unknown.The mutation types of the components are ordered by vertex labels.
WARNING:
- All finite types can be detected,
- All affine types can be detected, EXCEPT affine type D (the algorithm is not yet implemented)
- All exceptional types can be detected.
- Might fail to work if it is used within different Sage processes simultaneously (that happened in the doctesting).
EXAMPLES:
finite types:
sage: S = ClusterSeed(['A',5]) sage: S._mutation_type = S._quiver._mutation_type = None sage: S.mutation_type() ['A', 5] sage: S = ClusterSeed([(0,1),(1,2),(2,3),(3,4)]) sage: S.mutation_type() ['A', 5]
affine types:
sage: S = ClusterSeed(['E',8,[1,1]]); S A seed for a cluster algebra of rank 10 of type ['E', 8, [1, 1]] sage: S._mutation_type = S._quiver._mutation_type = None; S A seed for a cluster algebra of rank 10 sage: S.mutation_type() # long time ['E', 8, [1, 1]]
the not yet working affine type D:
sage: S = ClusterSeed(['D',4,1]) sage: S._mutation_type = S._quiver._mutation_type = None sage: S.mutation_type() # todo: not implemented ['D', 4, 1]
the exceptional types:
sage: S = ClusterSeed(['X',6]) sage: S._mutation_type = S._quiver._mutation_type = None sage: S.mutation_type() # long time ['X', 6]
infinite types:
sage: S = ClusterSeed(['GR',[4,9]]) sage: S._mutation_type = S._quiver._mutation_type = None sage: S.mutation_type() 'undetermined infinite mutation type'
-
mutations
()¶ Returns the list of mutations
self
has undergone if they are being tracked.Examples:
sage: S = ClusterSeed(['A',3]) sage: S.mutations() [] sage: S.mutate([0,1,0,2]) sage: S.mutations() [0, 1, 0, 2] sage: S.track_mutations(False) sage: S.mutations() Traceback (most recent call last): ... ValueError: Not recording mutation sequence. Need to track mutations.
-
n
()¶ Returns the number of exchangeable variables of
self
.EXAMPLES:
sage: S = ClusterSeed(['A',3]) sage: S.n() 3
-
oriented_exchange_graph
()¶ Return the oriented exchange graph of
self
as a directed graph.The seed must be a cluster seed for a cluster algebra of finite type with principal coefficients (the corresponding quiver must have mutable vertices 0,1,...,n-1).
EXAMPLES:
sage: S = ClusterSeed(['A', 2]).principal_extension() sage: G = S.oriented_exchange_graph(); G Digraph on 5 vertices sage: G.out_degree_sequence() [2, 1, 1, 1, 0] sage: S = ClusterSeed(['B', 2]).principal_extension() sage: G = S.oriented_exchange_graph(); G Digraph on 6 vertices sage: G.out_degree_sequence() [2, 1, 1, 1, 1, 0]
TESTS:
sage: S = ClusterSeed(['A',[2,2],1]) sage: S.oriented_exchange_graph() Traceback (most recent call last): ... TypeError: only works for finite mutation type sage: S = ClusterSeed(['A', 2]) sage: S.oriented_exchange_graph() Traceback (most recent call last): ... TypeError: only works for principal coefficients
-
plot
(circular=False, mark=None, save_pos=False, force_c=False, with_greens=False, add_labels=False)¶ Returns the plot of the quiver of
self
.INPUT:
circular
– (default:False) if True, the circular plot is chosen, otherwise >>spring<< is used.mark
– (default: None) if set to i, the vertex i is highlighted.save_pos
– (default:False) if True, the positions of the vertices are saved.force_c
– (default:False) if True, will show the frozen vertices even if they were never initializedwith_greens
– (default:False) if True, will display the green vertices in greenadd_labels
– (default:False) if True, will use the initial variables as labels
EXAMPLES:
sage: S = ClusterSeed(['A',5]) sage: pl = S.plot() sage: pl = S.plot(circular=True)
-
principal_extension
()¶ Returns the principal extension of self, yielding a 2n-by-n matrix. Raises an error if the input seed has a non-square exchange matrix. In this case, the method instead adds n frozen variables to any previously frozen variables. I.e., the seed obtained by adding a frozen variable to every exchangeable variable of
self
.EXAMPLES:
sage: S = ClusterSeed([[0,1],[1,2],[2,3],[2,4]]); S A seed for a cluster algebra of rank 5 sage: T = S.principal_extension(); T A seed for a cluster algebra of rank 5 with principal coefficients sage: T.b_matrix() [ 0 1 0 0 0] [-1 0 1 0 0] [ 0 -1 0 1 1] [ 0 0 -1 0 0] [ 0 0 -1 0 0] [ 1 0 0 0 0] [ 0 1 0 0 0] [ 0 0 1 0 0] [ 0 0 0 1 0] [ 0 0 0 0 1] sage: S = ClusterSeed(['A',4],user_labels=['a','b','c','d']) sage: T= S.principal_extension() sage: T.cluster() [a, b, c, d] sage: T.coefficients() [y0, y1, y2, y3] sage: S2 = ClusterSeed(['A',4],user_labels={0:'a',1:'b',2:'c',3:'d'}) sage: S2 == S True sage: T2 = S2.principal_extension() sage: T2 == T True
-
quiver
()¶ Returns the quiver associated to
self
.EXAMPLES:
sage: S = ClusterSeed(['A',3]) sage: S.quiver() Quiver on 3 vertices of type ['A', 3]
-
red_vertices
()¶ Return the list of red vertices of
self
.A vertex is defined to be red if its c-vector has all non-negative entries. More information on red vertices can be found at [BDP2013].
OUTPUT:
The red vertices as a list of integers.
EXAMPLES:
sage: ClusterSeed(['A',3]).principal_extension().red_vertices() [] sage: ClusterSeed(['A',[3,3],1]).principal_extension().red_vertices() [] sage: Q = ClusterSeed(['A',[3,3],1]).principal_extension(); sage: Q.mutate(1); sage: Q.red_vertices() [1]
-
reorient
(data)¶ Reorients
self
with respect to the given total order, or with respect to an iterator of ordered pairs.WARNING:
- This operation might change the mutation type of
self
. - Ignores ordered pairs \((i,j)\) for which neither \((i,j)\) nor \((j,i)\) is an edge of
self
.
INPUT:
data
– an iterator defining a total order onself.vertices()
, or an iterator of ordered pairs inself
defining the new orientation of these edges.
EXAMPLES:
sage: S = ClusterSeed(['A',[2,3],1]) sage: S.mutation_type() ['A', [2, 3], 1] sage: S.reorient([(0,1),(2,3)]) sage: S.mutation_type() ['D', 5] sage: S.reorient([(1,0),(2,3)]) sage: S.mutation_type() ['A', [1, 4], 1] sage: S.reorient([0,1,2,3,4]) sage: S.mutation_type() ['A', [1, 4], 1]
- This operation might change the mutation type of
-
reset_cluster
()¶ Resets the cluster of
self
to the initial cluster.EXAMPLES:
sage: S = ClusterSeed(['A',3]) sage: S.mutate([1,2,1]) sage: S.cluster() [x0, (x1 + 1)/x2, (x0*x2 + x1 + 1)/(x1*x2)] sage: S.reset_cluster() sage: S.cluster() [x0, x1, x2] sage: T = S.principal_extension() sage: T.cluster() [x0, x1, x2] sage: T.mutate([1,2,1]) sage: T.cluster() [x0, (x1*y2 + x0)/x2, (x1*y1*y2 + x0*y1 + x2)/(x1*x2)] sage: T.reset_cluster() sage: T.cluster() [x0, x1, x2] sage: S = ClusterSeed(['B',3],user_labels=[[1,2],[2,3],[3,4]],user_labels_prefix='p') sage: S.mutate([0,1]) sage: S.cluster() [(p_2_3 + 1)/p_1_2, (p_1_2*p_3_4^2 + p_2_3 + 1)/(p_1_2*p_2_3), p_3_4] sage: S.reset_cluster() sage: S.cluster() [p_1_2, p_2_3, p_3_4] sage: S.g_matrix() [1 0 0] [0 1 0] [0 0 1] sage: S.f_polynomials() [1, 1, 1]
-
reset_coefficients
()¶ Resets the coefficients of
self
to the frozen variables but keeps the current cluster. Raises an error if the number of frozen variables is different than the number of exchangeable variables.WARNING: This command to be phased out since ‘use_c_vectors() does this more effectively.
EXAMPLES:
sage: S = ClusterSeed(['A',3]).principal_extension() sage: S.b_matrix() [ 0 1 0] [-1 0 -1] [ 0 1 0] [ 1 0 0] [ 0 1 0] [ 0 0 1] sage: S.mutate([1,2,1]) sage: S.b_matrix() [ 0 1 -1] [-1 0 1] [ 1 -1 0] [ 1 0 0] [ 0 1 -1] [ 0 0 -1] sage: S.reset_coefficients() sage: S.b_matrix() [ 0 1 -1] [-1 0 1] [ 1 -1 0] [ 1 0 0] [ 0 1 0] [ 0 0 1]
-
save_image
(filename, circular=False, mark=None, save_pos=False)¶ Saves the plot of the underlying digraph of the quiver of
self
.INPUT:
filename
– the filename the image is saved to.circular
– (default: False) if True, the circular plot is chosen, otherwise >>spring<< is used.mark
– (default: None) if set to i, the vertex i is highlighted.save_pos
– (default:False) if True, the positions of the vertices are saved.
EXAMPLES:
sage: S = ClusterSeed(['F',4,[1,2]]) sage: S.save_image(os.path.join(SAGE_TMP, 'sage.png'))
-
set_c_matrix
(data)¶ Will force set the c matrix according to a matrix, a quiver, or a seed.
INPUT:
data
– The matrix to set the c matrix to. Also allowed to be a quiver or cluster seed, in which case the b_matrix is used.
EXAMPLES:
sage: S = ClusterSeed(['A',3]); sage: X = matrix([[0,0,1],[0,1,0],[1,0,0]]) sage: S.set_c_matrix(X) sage: S.c_matrix() [0 0 1] [0 1 0] [1 0 0] sage: Y = matrix([[-1,0,1],[0,1,0],[1,0,0]]) sage: S.set_c_matrix(Y) C matrix does not look to be valid - there exists a column containing positive and negative entries. Continuing... sage: Z = matrix([[1,0,1],[0,1,0],[2,0,2]]) sage: S.set_c_matrix(Z) C matrix does not look to be valid - not a linearly independent set. Continuing...
-
set_cluster
(cluster, force=False)¶ Sets the cluster for
self
tocluster
.Warning: Initialization may lead to inconsistent data.
INPUT:
cluster
– an iterable defining a cluster forself
.
EXAMPLES:
sage: S = ClusterSeed(['A',3]) sage: cluster = S.cluster() sage: S.mutate([1,2,1]) sage: S.cluster() [x0, (x1 + 1)/x2, (x0*x2 + x1 + 1)/(x1*x2)] sage: cluster2 = S.cluster() sage: S.set_cluster(cluster) Warning: using set_cluster at this point could lead to inconsistent seed data. sage: S.set_cluster(cluster, force=True) sage: S.cluster() [x0, x1, x2] sage: S.set_cluster(cluster2, force=True) sage: S.cluster() [x0, (x1 + 1)/x2, (x0*x2 + x1 + 1)/(x1*x2)] sage: S = ClusterSeed(['A',3]); S.use_fpolys(False) sage: S.set_cluster([1,1,1]) Warning: clusters not being tracked so this command is ignored.
-
show
(fig_size=1, circular=False, mark=None, save_pos=False, force_c=False, with_greens=False, add_labels=False)¶ Shows the plot of the quiver of
self
.INPUT:
fig_size
– (default: 1) factor by which the size of the plot is multiplied.circular
– (default: False) if True, the circular plot is chosen, otherwise >>spring<< is used.mark
– (default: None) if set to i, the vertex i is highlighted.save_pos
– (default:False) if True, the positions of the vertices are saved.force_c
– (default:False) if True, will show the frozen vertices even if they were never initializedwith_greens
– (default:False) if True, will display the green vertices in greenadd_labels
– (default:False) if True, will use the initial variables as labels
TESTS:
sage: S = ClusterSeed(['A',5]) sage: S.show() # long time
-
smallest_c_vector
()¶ Return the vertex with the smallest c vector
OUTPUT:
Returns an integer.
- EXAMPLES::
- sage: B = matrix([[0,2],[-2,0]]) sage: C = ClusterSeed(B).principal_extension(); sage: C.mutate(0) sage: C.smallest_c_vector() 0
-
track_mutations
(use=True)¶ Begins tracking the mutation path.
Warning: May initialize all other data to ensure that all c, d, and g vectors agree on the start of mutations.
INPUT:
use
– (default:True) If True, will begin filling the mutation path
EXAMPLES:
sage: S = ClusterSeed(['A',4]); S.track_mutations(False) sage: S.mutate(0) sage: S.mutations() Traceback (most recent call last): ... ValueError: Not recording mutation sequence. Need to track mutations. sage: S.track_mutations(True) sage: S.g_matrix() [1 0 0 0] [0 1 0 0] [0 0 1 0] [0 0 0 1] sage: S.mutate([0,1]) sage: S.mutations() [0, 1]
-
universal_extension
()¶ Returns the universal extension of
self
.This is the initial seed of the associated cluster algebra with universal coefficients, as defined in section 12 of Arxiv math/0602259.
This method works only if
self
is a bipartite, finite-type seed.Due to some limitations in the current implementation of
CartanType
, we need to construct the set of almost positive coroots by hand. As a consequence their ordering is not the standard one (the rows of the bottom part of the exchange matrix might be a shuffling of those you would expect).EXAMPLES:
sage: S = ClusterSeed(['A',2]) sage: T = S.universal_extension() sage: T.b_matrix() [ 0 1] [-1 0] [-1 0] [ 1 0] [ 1 -1] [ 0 1] [ 0 -1] sage: S = ClusterSeed(['A',3]) sage: T = S.universal_extension() sage: T.b_matrix() [ 0 1 0] [-1 0 -1] [ 0 1 0] [-1 0 0] [ 1 0 0] [ 1 -1 0] [ 1 -1 1] [ 0 1 0] [ 0 -1 0] [ 0 -1 1] [ 0 0 -1] [ 0 0 1] sage: S = ClusterSeed(['B',2]) sage: T = S.universal_extension() sage: T.b_matrix() [ 0 1] [-2 0] [-1 0] [ 1 0] [ 1 -1] [ 2 -1] [ 0 1] [ 0 -1]
-
urban_renewals
(return_first=False)¶ Return the list of the urban renewal vertices of
self
.An urban renewal vertex is one in which there are two arrows pointing toward the vertex and two arrows pointing away.
INPUT:
return_first
– (default:False) if True, will return the first urban renewal
OUTPUT:
A list of vertices (as integers)
EXAMPLES:
sage: G = ClusterSeed(['GR',[4,9]]); G.urban_renewals() [5, 6]
-
use_c_vectors
(use=True, bot_is_c=False, force=False)¶ Reconstruct c vectors from other data or initialize if no usable data exists.
Warning: Initialization may lead to inconsistent data.
INPUT:
use
– (default:True) If True, will use c vectorsbot_is_c
– (default:False) If True and ClusterSeed self has self._m == self._n, then will assume bottom half of the extended exchange matrix is the c-matrix. If true, lets the ClusterSeed know c-vectors can be calculated.
EXAMPLES:
sage: S = ClusterSeed(['A',4]); sage: S.use_c_vectors(False); S.use_g_vectors(False); S.use_fpolys(False); S.track_mutations(False) sage: S.use_c_vectors(True) Warning: Initializing c-vectors at this point could lead to inconsistent seed data. sage: S.use_c_vectors(True, force=True) sage: S.c_matrix() [1 0 0 0] [0 1 0 0] [0 0 1 0] [0 0 0 1] sage: S = ClusterSeed(['A',4]); sage: S.use_c_vectors(False); S.use_g_vectors(False); S.use_fpolys(False); S.track_mutations(False) sage: S.mutate(1); sage: S.use_c_vectors(True, force=True) sage: S.c_matrix() [1 0 0 0] [0 1 0 0] [0 0 1 0] [0 0 0 1]
-
use_d_vectors
(use=True, force=False)¶ Reconstruct d vectors from other data or initialize if no usable data exists.
Warning: Initialization may lead to inconsistent data.
INPUT:
use
– (default:True) If True, will use d vectors
EXAMPLES:
sage: S = ClusterSeed(['A',4]); sage: S.use_d_vectors(True) sage: S.d_matrix() [-1 0 0 0] [ 0 -1 0 0] [ 0 0 -1 0] [ 0 0 0 -1] sage: S = ClusterSeed(['A',4]); S.use_d_vectors(False); S.track_mutations(False); S.mutate(1); S.d_matrix() [-1 0 0 0] [ 0 1 0 0] [ 0 0 -1 0] [ 0 0 0 -1] sage: S.use_fpolys(False) sage: S.d_matrix() Traceback (most recent call last): ... ValueError: Unable to calculate d-vectors. Need to use d vectors. sage: S = ClusterSeed(['A',4]); S.use_d_vectors(False); S.track_mutations(False); S.mutate(1); S.d_matrix() [-1 0 0 0] [ 0 1 0 0] [ 0 0 -1 0] [ 0 0 0 -1] sage: S.use_fpolys(False) sage: S.use_d_vectors(True) Warning: Initializing d-vectors at this point could lead to inconsistent seed data. sage: S.use_d_vectors(True, force=True) sage: S.d_matrix() [-1 0 0 0] [ 0 -1 0 0] [ 0 0 -1 0] [ 0 0 0 -1] sage: S = ClusterSeed(['A',4]); S.mutate(1); S.d_matrix() [-1 0 0 0] [ 0 1 0 0] [ 0 0 -1 0] [ 0 0 0 -1] sage: S = ClusterSeed(['A',4]); S.use_d_vectors(True); S.mutate(1); S.d_matrix() [-1 0 0 0] [ 0 1 0 0] [ 0 0 -1 0] [ 0 0 0 -1]
-
use_fpolys
(use=True, user_labels=None, user_labels_prefix=None)¶ Use F-polynomials in our Cluster Seed
Note: This will automatically try to recompute the cluster variables if possible
INPUT:
use
– (default:True) If True, will use F-polynomialsuser_labels
– (default:None) If set will overwrite the default cluster variable labelsuser_labels_prefix
– (default:None) If set will overwrite the default
EXAMPLES:
sage: S = ClusterSeed(['A',4]); S.use_fpolys(False); S._cluster sage: S.use_fpolys(True) sage: S.cluster() [x0, x1, x2, x3] sage: S = ClusterSeed(['A',4]); S.use_fpolys(False); S.track_mutations(False); S.mutate(1) sage: S.use_fpolys(True) Traceback (most recent call last): ... ValueError: F-polynomials and Cluster Variables cannot be reconstructed from given data. sage: S.cluster() Traceback (most recent call last): ... ValueError: Clusters not being tracked
-
use_g_vectors
(use=True, force=False)¶ Reconstruct g vectors from other data or initialize if no usable data exists.
Warning: Initialization may lead to inconsistent data.
INPUT:
use
– (default:True) If True, will use g vectors
EXAMPLES:
sage: S = ClusterSeed(['A',4]); sage: S.use_g_vectors(False); S.use_fpolys(False) sage: S.use_g_vectors(True) sage: S.g_matrix() [1 0 0 0] [0 1 0 0] [0 0 1 0] [0 0 0 1] sage: S = ClusterSeed(['A',4]); sage: S.use_g_vectors(False); S.use_fpolys(False) sage: S.mutate(1); sage: S.use_g_vectors(True) sage: S.g_matrix() [ 1 0 0 0] [ 0 -1 0 0] [ 0 0 1 0] [ 0 0 0 1] sage: S = ClusterSeed(['A',4]); sage: S.use_g_vectors(False); S.use_fpolys(False); S.track_mutations(False) sage: S.mutate(1); sage: S.use_c_vectors(False) sage: S.g_matrix() Traceback (most recent call last): ... ValueError: Unable to calculate g-vectors. Need to use g vectors. sage: S = ClusterSeed(['A',4]); sage: S.use_g_vectors(False); S.use_fpolys(False); S.track_mutations(False) sage: S.mutate(1); sage: S.use_c_vectors(False) sage: S.use_g_vectors(True) Warning: Initializing g-vectors at this point could lead to inconsistent seed data. sage: S.use_g_vectors(True, force=True) sage: S.g_matrix() [1 0 0 0] [0 1 0 0] [0 0 1 0] [0 0 0 1]
-
variable_class
(depth=+Infinity, ignore_bipartite_belt=False)¶ Returns all cluster variables in the mutation class of
self
.INPUT:
depth
– (default:infinity) integer, only seeds with distance at most depth from self are returnedignore_bipartite_belt
– (default:False) if True, the algorithms does not use the bipartite belt
EXAMPLES:
- for examples see
variable_class_iter()
TESTS:
sage: A = ClusterSeed(['A',3]).variable_class()
-
variable_class_iter
(depth=+Infinity, ignore_bipartite_belt=False)¶ Returns an iterator for all cluster variables in the mutation class of
self
.INPUT:
depth
– (default:infinity) integer, only seeds with distance at most depth from self are returnedignore_bipartite_belt
– (default:False) if True, the algorithms does not use the bipartite belt
EXAMPLES:
A standard finite type example:
sage: S = ClusterSeed(['A',3]) sage: it = S.variable_class_iter() sage: for T in it: print(T) x0 x1 x2 (x1 + 1)/x0 (x1^2 + x0*x2 + 2*x1 + 1)/(x0*x1*x2) (x1 + 1)/x2 (x0*x2 + x1 + 1)/(x0*x1) (x0*x2 + 1)/x1 (x0*x2 + x1 + 1)/(x1*x2)
Finite type examples with given depth:
sage: it = S.variable_class_iter(depth=1) sage: for T in it: print(T) Found a bipartite seed - restarting the depth counter at zero and constructing the variable class using its bipartite belt. x0 x1 x2 (x1 + 1)/x0 (x1^2 + x0*x2 + 2*x1 + 1)/(x0*x1*x2) (x1 + 1)/x2 (x0*x2 + x1 + 1)/(x0*x1) (x0*x2 + 1)/x1 (x0*x2 + x1 + 1)/(x1*x2)
Note that the notion of depth depends on whether a bipartite seed is found or not, or if it is manually ignored:
sage: it = S.variable_class_iter(depth=1,ignore_bipartite_belt=True) sage: for T in it: print(T) x0 x1 x2 (x1 + 1)/x2 (x0*x2 + 1)/x1 (x1 + 1)/x0 sage: S.mutate([0,1]) sage: it2 = S.variable_class_iter(depth=1) sage: for T in it2: print(T) (x1 + 1)/x0 (x0*x2 + x1 + 1)/(x0*x1) x2 (x1^2 + x0*x2 + 2*x1 + 1)/(x0*x1*x2) x1 (x0*x2 + 1)/x1
Infinite type examples:
sage: S = ClusterSeed(['A',[1,1],1]) sage: it = S.variable_class_iter(depth=2) sage: for T in it: print(T) Found a bipartite seed - restarting the depth counter at zero and constructing the variable class using its bipartite belt. x0 x1 (x1^2 + 1)/x0 (x1^4 + x0^2 + 2*x1^2 + 1)/(x0^2*x1) (x0^4 + 2*x0^2 + x1^2 + 1)/(x0*x1^2) (x0^2 + 1)/x1 (x1^6 + x0^4 + 2*x0^2*x1^2 + 3*x1^4 + 2*x0^2 + 3*x1^2 + 1)/(x0^3*x1^2) (x1^8 + x0^6 + 2*x0^4*x1^2 + 3*x0^2*x1^4 + 4*x1^6 + 3*x0^4 + 6*x0^2*x1^2 + 6*x1^4 + 3*x0^2 + 4*x1^2 + 1)/(x0^4*x1^3) (x0^8 + 4*x0^6 + 3*x0^4*x1^2 + 2*x0^2*x1^4 + x1^6 + 6*x0^4 + 6*x0^2*x1^2 + 3*x1^4 + 4*x0^2 + 3*x1^2 + 1)/(x0^3*x1^4) (x0^6 + 3*x0^4 + 2*x0^2*x1^2 + x1^4 + 3*x0^2 + 2*x1^2 + 1)/(x0^2*x1^3)
-
x
(k)¶ Returns the \(k\) -th initial cluster variable for the associated cluster seed.
EXAMPLES:
sage: S = ClusterSeed(['A',3]) sage: S.mutate([2,1]) sage: S.x(0) x0 sage: S.x(1) x1 sage: S.x(2) x2
-
y
(k)¶ Returns the \(k\) -th initial coefficient (frozen variable) for the associated cluster seed.
EXAMPLES:
sage: S = ClusterSeed(['A',3]).principal_extension() sage: S.mutate([2,1]) sage: S.y(0) y0 sage: S.y(1) y1 sage: S.y(2) y2
-
class
sage.combinat.cluster_algebra_quiver.cluster_seed.
ClusterVariable
(parent, numerator, denominator, coerce=True, reduce=True, mutation_type=None, variable_type=None, xdim=0)¶ Bases:
sage.rings.fraction_field_element.FractionFieldElement
This class is a thin wrapper for cluster variables in cluster seeds.
It provides the extra feature to store if a variable is frozen or not.
the associated positive root:
sage: S = ClusterSeed(['A',3]) sage: for T in S.variable_class_iter(): ....: print("{} {}".format(T, T.almost_positive_root())) x0 -alpha[1] x1 -alpha[2] x2 -alpha[3] (x1 + 1)/x0 alpha[1] (x1^2 + x0*x2 + 2*x1 + 1)/(x0*x1*x2) alpha[1] + alpha[2] + alpha[3] (x1 + 1)/x2 alpha[3] (x0*x2 + x1 + 1)/(x0*x1) alpha[1] + alpha[2] (x0*x2 + 1)/x1 alpha[2] (x0*x2 + x1 + 1)/(x1*x2) alpha[2] + alpha[3]
-
almost_positive_root
()¶ Returns the almost positive root associated to
self
ifself
is of finite type.EXAMPLES:
sage: S = ClusterSeed(['A',3]) sage: for T in S.variable_class_iter(): ....: print("{} {}".format(T, T.almost_positive_root())) x0 -alpha[1] x1 -alpha[2] x2 -alpha[3] (x1 + 1)/x0 alpha[1] (x1^2 + x0*x2 + 2*x1 + 1)/(x0*x1*x2) alpha[1] + alpha[2] + alpha[3] (x1 + 1)/x2 alpha[3] (x0*x2 + x1 + 1)/(x0*x1) alpha[1] + alpha[2] (x0*x2 + 1)/x1 alpha[2] (x0*x2 + x1 + 1)/(x1*x2) alpha[2] + alpha[3]
-
sage.combinat.cluster_algebra_quiver.cluster_seed.
PathSubset
(n, m)¶ Encodes a maximal Dyck path from (0,0) to (n,m) (for n >= m >= 0) as a subset of {0,1,2,..., 2n-1}. The encoding is given by indexing horizontal edges by odd numbers and vertical edges by evens.
The horizontal between (i,j) and (i+1,j) is indexed by the odd number 2*i+1. The vertical between (i,j) and (i,j+1) is indexed by the even number 2*j.
EXAMPLES:
sage: from sage.combinat.cluster_algebra_quiver.cluster_seed import PathSubset sage: PathSubset(4,0) {1, 3, 5, 7} sage: PathSubset(4,1) {1, 3, 5, 6, 7} sage: PathSubset(4,2) {1, 2, 3, 5, 6, 7} sage: PathSubset(4,3) {1, 2, 3, 4, 5, 6, 7} sage: PathSubset(4,4) {0, 1, 2, 3, 4, 5, 6, 7}
-
sage.combinat.cluster_algebra_quiver.cluster_seed.
SetToPath
(T)¶ Rearranges the encoding for a maximal Dyck path (as a set) so that it is a list in the proper order of the edges.
EXAMPLES:
sage: from sage.combinat.cluster_algebra_quiver.cluster_seed import PathSubset sage: from sage.combinat.cluster_algebra_quiver.cluster_seed import SetToPath sage: SetToPath(PathSubset(4,0)) [1, 3, 5, 7] sage: SetToPath(PathSubset(4,1)) [1, 3, 5, 7, 6] sage: SetToPath(PathSubset(4,2)) [1, 3, 2, 5, 7, 6] sage: SetToPath(PathSubset(4,3)) [1, 3, 2, 5, 4, 7, 6] sage: SetToPath(PathSubset(4,4)) [1, 0, 3, 2, 5, 4, 7, 6]
-
sage.combinat.cluster_algebra_quiver.cluster_seed.
coeff_recurs
(p, q, a1, a2, b, c)¶ Coefficients in Laurent expansion of greedy element, as defined by recursion.
EXAMPLES:
sage: from sage.combinat.cluster_algebra_quiver.cluster_seed import coeff_recurs sage: coeff_recurs(1, 1, 5, 5, 3, 3) 10
-
sage.combinat.cluster_algebra_quiver.cluster_seed.
get_green_vertices
(C)¶ Get the green vertices from a matrix. Will go through each clumn and return the ones where no entry is greater than 0.
INPUT:
C
– The C matrix to check
EXAMPLES:
sage: from sage.combinat.cluster_algebra_quiver.cluster_seed import get_green_vertices sage: S = ClusterSeed(['A',4]); S.mutate([1,2,3,2,0,1,2,0,3]) sage: get_green_vertices(S.c_matrix()) [0, 3]
-
sage.combinat.cluster_algebra_quiver.cluster_seed.
get_red_vertices
(C)¶ Get the red vertices from a matrix. Will go through each clumn and return the ones where no entry is less than 0.
INPUT:
C
– The C matrix to check
- EXAMPLES::
- sage: from sage.combinat.cluster_algebra_quiver.cluster_seed import get_red_vertices sage: S = ClusterSeed([‘A’,4]); S.mutate([1,2,3,2,0,1,2,0,3]) sage: get_red_vertices(S.c_matrix()) [1, 2]
-
sage.combinat.cluster_algebra_quiver.cluster_seed.
is_LeeLiZel_allowable
(T, n, m, b, c)¶ Check if the subset T contributes to the computation of the greedy element x[m,n] in the rank two (b,c)-cluster algebra.
This uses the conditions of Lee-Li-Zelevinsky’s paper [LeeLiZe].
EXAMPLES:
sage: from sage.combinat.cluster_algebra_quiver.cluster_seed import is_LeeLiZel_allowable sage: is_LeeLiZel_allowable({1,3,2,5,7,6},4,2,6,6) False sage: is_LeeLiZel_allowable({1,2,5},3,3,1,1) True