Hot-keys on this page
r m x p toggle line displays
j k next/prev highlighted chunk
0 (zero) top of page
1 (one) first highlighted chunk
r""" Groebner Fans
Sage provides much of the functionality of gfan, which is a software package whose main function is to enumerate all reduced Groebner bases of a polynomial ideal. The reduced Groebner bases yield the maximal cones in the Groebner fan of the ideal. Several subcomputations can be issued and additional tools are included. Among these the highlights are:
- Commands for computing tropical varieties.
- Interactive walks in the Groebner fan of an ideal.
- Commands for graphical renderings of Groebner fans and monomial ideals.
AUTHORS:
- Anders Nedergaard Jensen: Wrote the gfan C++ program, which implements algorithms many of which were invented by Jensen, Komei Fukuda, and Rekha Thomas. All the underlying hard work of the Groebner fans functionality of Sage depends on this C++ program.
- William Stein (2006-04-20): Wrote first version of the Sage code for working with Groebner fans.
- Tristram Bogart: the design of the Sage interface to gfan is joint work with Tristram Bogart, who also supplied numerous examples.
- Marshall Hampton (2008-03-25): Rewrote various functions to use gfan-0.3. This is still a work in progress, comments are appreciated on sage-devel@googlegroups.com (or personally at hamptonio@gmail.com).
EXAMPLES::
sage: x,y = QQ['x,y'].gens() sage: i = ideal(x^2 - y^2 + 1) sage: g = i.groebner_fan() sage: g.reduced_groebner_bases() [[x^2 - y^2 + 1], [-x^2 + y^2 - 1]]
TESTS::
sage: x,y = QQ['x,y'].gens() sage: i = ideal(x^2 - y^2 + 1) sage: g = i.groebner_fan() sage: g == loads(dumps(g)) True
REFERENCES:
- Anders N. Jensen; Gfan, a software system for Groebner fans; available at http://www.math.tu-berlin.de/~jensen/software/gfan/gfan.html """
""" Checks if any strings in a list are prefixes of another string in the list.
EXAMPLES::
sage: from sage.rings.polynomial.groebner_fan import prefix_check sage: prefix_check(['z1','z1z1']) False sage: prefix_check(['z1','zz1']) True """
""" Computes the maximum degree of a list of polynomials
EXAMPLES::
sage: from sage.rings.polynomial.groebner_fan import max_degree sage: R.<x,y> = PolynomialRing(QQ,2) sage: p_list = [x^2-y,x*y^10-x] sage: max_degree(p_list) 11.0 """
""" Utility function that parses cone information into a dict indexed by dimension.
INPUT:
- ``fan_dict_cone`` - the value of a fan_dict with key 'CONES'
EXAMPLES::
sage: R.<x,y,z,w> = PolynomialRing(QQ,4) sage: I = R.ideal([x^2+y^2+z^2-1,x^4-y^2-z-1,x+y+z,w+x+y]) sage: GF = I.groebner_fan() sage: TI = GF.tropical_intersection() sage: from sage.rings.polynomial.groebner_fan import _cone_parse sage: _cone_parse(TI.fan_dict['CONES']) {1: [[0], [1], [2], [3], [4]], 2: [[2, 3]]} """ temp = [''] else: else:
""" Converts polymake/gfan data on a polyhedral cone into a sage class. Currently (18-03-2008) needs a lot of work.
EXAMPLES::
sage: R3.<x,y,z> = PolynomialRing(QQ,3) sage: gf = R3.ideal([x^8-y^4,y^4-z^2,z^2-2]).groebner_fan() sage: a = gf[0].groebner_cone() sage: a.facets() [[0, 0, 1], [0, 1, 0], [1, 0, 0]] """
""" Returns a basic description of the polyhedral cone.
EXAMPLES::
sage: R3.<x,y,z> = PolynomialRing(QQ,3) sage: gf = R3.ideal([x^8-y^4,y^4-z^2,z^2-2]).groebner_fan() sage: a = gf[0].groebner_cone() sage: a # indirect doctests Polyhedral cone in 3 dimensions of dimension 3 """
""" Returns the inward facet normals of the Groebner cone.
EXAMPLES::
sage: R3.<x,y,z> = PolynomialRing(QQ,3) sage: gf = R3.ideal([x^8-y^4,y^4-z^2,z^2-2]).groebner_fan() sage: a = gf[0].groebner_cone() sage: a.facets() [[0, 0, 1], [0, 1, 0], [1, 0, 0]] """
""" Returns the ambient dimension of the Groebner cone.
EXAMPLES::
sage: R3.<x,y,z> = PolynomialRing(QQ,3) sage: gf = R3.ideal([x^8-y^4,y^4-z^2,z^2-2]).groebner_fan() sage: a = gf[0].groebner_cone() sage: a.ambient_dim() 3 """
""" Returns the dimension of the Groebner cone.
EXAMPLES::
sage: R3.<x,y,z> = PolynomialRing(QQ,3) sage: gf = R3.ideal([x^8-y^4,y^4-z^2,z^2-2]).groebner_fan() sage: a = gf[0].groebner_cone() sage: a.dim() 3 """
""" Returns the lineality dimension of the Groebner cone. This is just the difference between the ambient dimension and the dimension of the cone.
EXAMPLES::
sage: R3.<x,y,z> = PolynomialRing(QQ,3) sage: gf = R3.ideal([x^8-y^4,y^4-z^2,z^2-2]).groebner_fan() sage: a = gf[0].groebner_cone() sage: a.lineality_dim() 0 """
""" Returns a point in the relative interior of the Groebner cone.
EXAMPLES::
sage: R3.<x,y,z> = PolynomialRing(QQ,3) sage: gf = R3.ideal([x^8-y^4,y^4-z^2,z^2-2]).groebner_fan() sage: a = gf[0].groebner_cone() sage: a.relative_interior_point() [1, 1, 1] """
""" Converts polymake/gfan data on a polyhedral fan into a sage class.
INPUT:
- ``gfan_polyhedral_fan`` - output from gfan of a polyhedral fan.
EXAMPLES::
sage: R.<x,y,z> = PolynomialRing(QQ,3) sage: i2 = ideal(x*z + 6*y*z - z^2, x*y + 6*x*z + y*z - z^2, y^2 + x*z + y*z) sage: gf2 = i2.groebner_fan(verbose = False) sage: pf = gf2.polyhedralfan() sage: pf.rays() [[-1, 0, 1], [-1, 1, 0], [1, -2, 1], [1, 1, -2], [2, -1, -1]] """ 'LINEALITY_SPACE','ORTH_LINEALITY_SPACE','F_VECTOR', 'CONES','MAXIMAL_CONES','PURE','SIMPLICIAL','MULTIPLICITIES']
""" Returns a basic description of the polyhedral fan.
EXAMPLES::
sage: R3.<x,y,z> = PolynomialRing(QQ,3) sage: gf = R3.ideal([x^8-y^4,y^4-z^2,z^2-2]).groebner_fan() sage: pf = gf.polyhedralfan() sage: pf # indirect doctest Polyhedral fan in 3 dimensions of dimension 3 """
r""" Returns the raw output of gfan as a string. This should only be needed internally as all relevant output is converted to sage objects.
EXAMPLES::
sage: R3.<x,y,z> = PolynomialRing(QQ,3) sage: gf = R3.ideal([x^8-y^4,y^4-z^2,z^2-2]).groebner_fan() sage: pf = gf.polyhedralfan() sage: pf._str_() '_application fan\n_version 2.2\n_type SymmetricFan\n\nAMBIENT_DIM\n3\n\nDIM\n3\n\nLINEALITY_DIM\n0\n\nRAYS\n0 0 1\t# 0\n0 1 0\t# 1\n1 0 0\t# 2\n\nN_RAYS\n3\n\nLINEALITY_SPACE\n\nORTH_LINEALITY_SPACE\n1 0 0\n0 1 0\n0 0 1\n\nF_VECTOR\n1 3 3 1\n\nSIMPLICIAL\n1\n\nPURE\n1\n\nCONES\n{}\t# Dimension 0\n{0}\t# Dimension 1\n{1}\n{2}\n{0 1}\t# Dimension 2\n{0 2}\n{1 2}\n{0 1 2}\t# Dimension 3\n\nMAXIMAL_CONES\n{0 1 2}\t# Dimension 3\n' """
""" Returns the ambient dimension of the Groebner fan.
EXAMPLES::
sage: R3.<x,y,z> = PolynomialRing(QQ,3) sage: gf = R3.ideal([x^8-y^4,y^4-z^2,z^2-2]).groebner_fan() sage: a = gf.polyhedralfan() sage: a.ambient_dim() 3 """
""" Returns the dimension of the Groebner fan.
EXAMPLES::
sage: R3.<x,y,z> = PolynomialRing(QQ,3) sage: gf = R3.ideal([x^8-y^4,y^4-z^2,z^2-2]).groebner_fan() sage: a = gf.polyhedralfan() sage: a.dim() 3 """
""" Returns the lineality dimension of the fan. This is the dimension of the largest subspace contained in the fan.
EXAMPLES::
sage: R3.<x,y,z> = PolynomialRing(QQ,3) sage: gf = R3.ideal([x^8-y^4,y^4-z^2,z^2-2]).groebner_fan() sage: a = gf.polyhedralfan() sage: a.lineality_dim() 0 """
""" A list of rays of the polyhedral fan.
EXAMPLES::
sage: R.<x,y,z> = PolynomialRing(QQ,3) sage: i2 = ideal(x*z + 6*y*z - z^2, x*y + 6*x*z + y*z - z^2, y^2 + x*z + y*z) sage: gf2 = i2.groebner_fan(verbose = False) sage: pf = gf2.polyhedralfan() sage: pf.rays() [[-1, 0, 1], [-1, 1, 0], [1, -2, 1], [1, 1, -2], [2, -1, -1]] """
""" A dictionary of cones in which the keys are the cone dimensions. For each dimension, the value is a list of the cones, where each element consists of a list of ray indices.
EXAMPLES::
sage: R.<x,y,z> = QQ[] sage: f = 1+x+y+x*y sage: I = R.ideal([f+z*f, 2*f+z*f, 3*f+z^2*f]) sage: GF = I.groebner_fan() sage: PF = GF.tropical_intersection() sage: PF.cones() {1: [[0], [1], [2], [3], [4], [5]], 2: [[0, 1], [0, 2], [0, 3], [0, 4], [1, 2], [1, 3], [2, 4], [3, 4], [1, 5], [2, 5], [3, 5], [4, 5]]} """
""" A dictionary of the maximal cones in which the keys are the cone dimensions. For each dimension, the value is a list of the maximal cones, where each element consists of a list of ray indices.
EXAMPLES::
sage: R.<x,y,z> = QQ[] sage: f = 1+x+y+x*y sage: I = R.ideal([f+z*f, 2*f+z*f, 3*f+z^2*f]) sage: GF = I.groebner_fan() sage: PF = GF.tropical_intersection() sage: PF.maximal_cones() {2: [[0, 1], [0, 2], [0, 3], [0, 4], [1, 2], [1, 3], [2, 4], [3, 4], [1, 5], [2, 5], [3, 5], [4, 5]]} """
""" The f-vector of the fan.
EXAMPLES::
sage: R.<x,y,z> = QQ[] sage: f = 1+x+y+x*y sage: I = R.ideal([f+z*f, 2*f+z*f, 3*f+z^2*f]) sage: GF = I.groebner_fan() sage: PF = GF.tropical_intersection() sage: PF.f_vector() [1, 6, 12] """
""" Whether the fan is simplicial or not.
EXAMPLES::
sage: R.<x,y,z> = QQ[] sage: f = 1+x+y+x*y sage: I = R.ideal([f+z*f, 2*f+z*f, 3*f+z^2*f]) sage: GF = I.groebner_fan() sage: PF = GF.tropical_intersection() sage: PF.is_simplicial() True """
""" Converts to the RationalPolyhedralFan class, which is more actively maintained. While the information in each class is essentially the same, the methods and implementation are different.
EXAMPLES::
sage: R.<x,y,z> = QQ[] sage: f = 1+x+y+x*y sage: I = R.ideal([f+z*f, 2*f+z*f, 3*f+z^2*f]) sage: GF = I.groebner_fan() sage: PF = GF.tropical_intersection() sage: fan = PF.to_RationalPolyhedralFan() sage: [tuple(q.facet_normals()) for q in fan] [(M(0, -1, 0), M(-1, 0, 0)), (M(0, 0, -1), M(-1, 0, 0)), (M(0, 0, 1), M(-1, 0, 0)), (M(0, 1, 0), M(-1, 0, 0)), (M(0, 0, -1), M(0, -1, 0)), (M(0, 0, 1), M(0, -1, 0)), (M(0, 1, 0), M(0, 0, -1)), (M(0, 1, 0), M(0, 0, 1)), (M(1, 0, 0), M(0, -1, 0)), (M(1, 0, 0), M(0, 0, -1)), (M(1, 0, 0), M(0, 0, 1)), (M(1, 0, 0), M(0, 1, 0))]
Here we use the RationalPolyhedralFan's Gale_transform method on a tropical prevariety.
.. link
::
sage: fan.Gale_transform() [ 1 0 0 0 0 1 -2] [ 0 1 0 0 1 0 -2] [ 0 0 1 1 0 0 -2]
"""
""" A system of initial forms from a polynomial system. To each form is associated a cone and a list of polynomials (the initial form system itself).
This class is intended for internal use inside of the TropicalPrevariety class.
EXAMPLES::
sage: from sage.rings.polynomial.groebner_fan import InitialForm sage: R.<x,y> = QQ[] sage: inform = InitialForm([0], [[-1, 0]], [y^2 - 1, y^2 - 2, y^2 - 3]) sage: inform._cone [0] """
""" The cone associated with the initial form system.
EXAMPLES::
sage: R.<x,y> = QQ[] sage: I = R.ideal([(x+y)^2-1,(x+y)^2-2,(x+y)^2-3]) sage: GF = I.groebner_fan() sage: PF = GF.tropical_intersection() sage: pfi0 = PF.initial_form_systems()[0] sage: pfi0.cone() [0] """
""" The rays of the cone associated with the initial form system.
EXAMPLES::
sage: R.<x,y> = QQ[] sage: I = R.ideal([(x+y)^2-1,(x+y)^2-2,(x+y)^2-3]) sage: GF = I.groebner_fan() sage: PF = GF.tropical_intersection() sage: pfi0 = PF.initial_form_systems()[0] sage: pfi0.rays() [[-1, 0]] """
""" A ray internal to the cone associated with the initial form system.
EXAMPLES::
sage: R.<x,y> = QQ[] sage: I = R.ideal([(x+y)^2-1,(x+y)^2-2,(x+y)^2-3]) sage: GF = I.groebner_fan() sage: PF = GF.tropical_intersection() sage: pfi0 = PF.initial_form_systems()[0] sage: pfi0.internal_ray() (-1, 0) """
""" The initial forms (polynomials).
EXAMPLES::
sage: R.<x,y> = QQ[] sage: I = R.ideal([(x+y)^2-1,(x+y)^2-2,(x+y)^2-3]) sage: GF = I.groebner_fan() sage: PF = GF.tropical_intersection() sage: pfi0 = PF.initial_form_systems()[0] sage: pfi0.initial_forms() [y^2 - 1, y^2 - 2, y^2 - 3] """
""" Returns the exponents of the vertices of a Newton polytope that make up the supporting hyperplane for the given outward normal.
EXAMPLES::
sage: from sage.rings.polynomial.groebner_fan import verts_for_normal sage: R.<x,y,z> = PolynomialRing(QQ,3) sage: f1 = x*y*z - 1 sage: f2 = f1*(x^2 + y^2 + 1) sage: verts_for_normal([1,1,1],f2) [(3, 1, 1), (1, 3, 1)] """
""" This class is a subclass of the PolyhedralFan class, with some additional methods for tropical prevarieties.
INPUT:
- ``gfan_polyhedral_fan`` - output from gfan of a polyhedral fan. - ``polynomial_system`` - a list of polynomials - ``poly_ring`` - the polynomial ring of the list of polynomials - ``parameters`` (optional) - a list of variables to be considered as parameters
EXAMPLES::
sage: R.<x,y,z> = QQ[] sage: I = R.ideal([(x+y+z)^2-1,(x+y+z)-x,(x+y+z)-3]) sage: GF = I.groebner_fan() sage: TI = GF.tropical_intersection() sage: TI._polynomial_system [x^2 + 2*x*y + y^2 + 2*x*z + 2*y*z + z^2 - 1, y + z, x + y + z - 3] """
""" Returns a list of systems of initial forms for each cone in the tropical prevariety.
EXAMPLES::
sage: R.<x,y> = QQ[] sage: I = R.ideal([(x+y)^2-1,(x+y)^2-2,(x+y)^2-3]) sage: GF = I.groebner_fan() sage: PF = GF.tropical_intersection() sage: pfi = PF.initial_form_systems() sage: for q in pfi: ....: print(q.initial_forms()) [y^2 - 1, y^2 - 2, y^2 - 3] [x^2 - 1, x^2 - 2, x^2 - 3] [x^2 + 2*x*y + y^2, x^2 + 2*x*y + y^2, x^2 + 2*x*y + y^2] """
""" Converts a ring to gfan's format.
EXAMPLES::
sage: R.<w,x,y,z> = QQ[] sage: from sage.rings.polynomial.groebner_fan import ring_to_gfan_format sage: ring_to_gfan_format(R) 'Q[w, x, y, z]' sage: R2.<x,y> = GF(2)[] sage: ring_to_gfan_format(R2) 'Z/2Z[x, y]' """ ring_str = 'Z' + str(input_ring.gens()).replace('(','[').replace(')',']') else:
""" Return the ideal in gfan's notation.
EXAMPLES::
sage: R.<x,y,z> = PolynomialRing(QQ,3) sage: polys = [x^2*y - z, y^2*z - x, z^2*x - y] sage: from sage.rings.polynomial.groebner_fan import ideal_to_gfan_format sage: ideal_to_gfan_format(R, polys) 'Q[x, y, z]{x^2*y-z,y^2*z-x,x*z^2-y}'
TESTS:
Test that :trac:`20146` is fixed::
sage: P = PolynomialRing(QQ,"x11,x12,x13,x14,x15,x21,x22,x23,x24,x25,x31,x32,x33,x34,x35"); x = P.gens(); M = Matrix(3,x) sage: I = P.ideal(M.minors(2)) sage: ideal_to_gfan_format(P,I.gens()) 'Q[x11, x12, x13, x14, x15, x21, x22, x23, x24, x25, x31, x32, x33, x34, x35]{-x12*x21+x11*x22,-x13*x21+x11*x23,-x14*x21+x11*x24,-x15*x21+x11*x25,-x13*x22+x12*x23,-x14*x22+x12*x24,-x15*x22+x12*x25,-x14*x23+x13*x24,-x15*x23+x13*x25,-x15*x24+x14*x25,-x12*x31+x11*x32,-x13*x31+x11*x33,-x14*x31+x11*x34,-x15*x31+x11*x35,-x13*x32+x12*x33,-x14*x32+x12*x34,-x15*x32+x12*x35,-x14*x33+x13*x34,-x15*x33+x13*x35,-x15*x34+x14*x35,-x22*x31+x21*x32,-x23*x31+x21*x33,-x24*x31+x21*x34,-x25*x31+x21*x35,-x23*x32+x22*x33,-x24*x32+x22*x34,-x25*x32+x22*x35,-x24*x33+x23*x34,-x25*x33+x23*x35,-x25*x34+x24*x35}' """
""" This class is used to access capabilities of the program Gfan. In addition to computing Groebner fans, Gfan can compute other things in tropical geometry such as tropical prevarieties.
INPUT:
- ``I`` - ideal in a multivariate polynomial ring
- ``is_groebner_basis`` - bool (default False). if True, then I.gens() must be a Groebner basis with respect to the standard degree lexicographic term order.
- ``symmetry`` - default: None; if not None, describes symmetries of the ideal
- ``verbose`` - default: False; if True, printout useful info during computations
EXAMPLES::
sage: R.<x,y,z> = QQ[] sage: I = R.ideal([x^2*y - z, y^2*z - x, z^2*x - y]) sage: G = I.groebner_fan(); G Groebner fan of the ideal: Ideal (x^2*y - z, y^2*z - x, x*z^2 - y) of Multivariate Polynomial Ring in x, y, z over Rational Field
Here is an example of the use of the tropical_intersection command, and then using the RationalPolyhedralFan class to compute the Stanley-Reisner ideal of the tropical prevariety::
sage: R.<x,y,z> = QQ[] sage: I = R.ideal([(x+y+z)^3-1,(x+y+z)^3-x,(x+y+z)-3]) sage: GF = I.groebner_fan() sage: PF = GF.tropical_intersection() sage: PF.rays() [[-1, 0, 0], [0, -1, 0], [0, 0, -1], [1, 1, 1]] sage: RPF = PF.to_RationalPolyhedralFan() sage: RPF.Stanley_Reisner_ideal(PolynomialRing(QQ,4,'A, B, C, D')) Ideal (A*B, A*C, B*C*D) of Multivariate Polynomial Ring in A, B, C, D over Rational Field
""" print("WARNING! Symmetry option not yet implemented!!") raise TypeError("I must be a multivariate polynomial ideal") raise RuntimeError("Ring variables cannot contain each other as prefixes") # todo: add support for ZZ, which only works for bases computation, not tropical intersections raise NotImplementedError("Groebner fan computation only implemented over fields") # sage: previous_prime (2^15) # 32749 raise NotImplementedError("Groebner fan computation only implemented over Q or GF(p) for p <= 32749.") raise NotImplementedError("Groebner fan computation only implemented for rings in at most 52 variables.")
""" Describes the Groebner fan and its corresponding ideal.
EXAMPLES::
sage: R.<q,u> = PolynomialRing(QQ,2) sage: gf = R.ideal([q-u,u^2-1]).groebner_fan() sage: gf # indirect doctest Groebner fan of the ideal: Ideal (q - u, u^2 - 1) of Multivariate Polynomial Ring in q, u over Rational Field """
""" Tests equality of Groebner fan objects.
EXAMPLES::
sage: R.<q,u> = PolynomialRing(QQ,2) sage: gf = R.ideal([q^2-u,u^2-q]).groebner_fan() sage: gf2 = R.ideal([u^2-q,q^2-u]).groebner_fan() sage: gf.__eq__(gf2) True """
""" Return the ideal the was used to define this Groebner fan.
EXAMPLES::
sage: R.<x1,x2> = PolynomialRing(QQ,2) sage: gf = R.ideal([x1^3-x2,x2^3-2*x1-2]).groebner_fan() sage: gf.ideal() Ideal (x1^3 - x2, x2^3 - 2*x1 - 2) of Multivariate Polynomial Ring in x1, x2 over Rational Field """
""" INPUT: none
OUTPUT:
- map from Sage ring to gfan ring
- map from gfan ring to Sage ring
EXAMPLES::
sage: R.<x,y,z> = PolynomialRing(QQ,3) sage: G = R.ideal([x^2*y - z, y^2*z - x, z^2*x - y]).groebner_fan() sage: G._gfan_maps() (Ring morphism: From: Multivariate Polynomial Ring in x, y, z over Rational Field To: Multivariate Polynomial Ring in a, b, c over Rational Field Defn: x |--> a y |--> b z |--> c, Ring morphism: From: Multivariate Polynomial Ring in a, b, c over Rational Field To: Multivariate Polynomial Ring in x, y, z over Rational Field Defn: a |--> x b |--> y c |--> z) """
# Define a polynomial ring in n variables # that are named a,b,c,d, ..., z, A, B, C, ...
# Define the homomorphism that sends the # generators of S to the generators of T.
# Define the homomorphism that sends the # generators of T to the generators of S.
""" Return the ring in gfan's notation
EXAMPLES::
sage: R.<x,y,z> = PolynomialRing(QQ,3) sage: G = R.ideal([x^2*y - z, y^2*z - x, z^2*x - y]).groebner_fan() sage: G._gfan_ring() 'Q[x, y, z]' """
""" Return the ideal in gfan's notation.
EXAMPLES::
sage: R.<x,y,z> = PolynomialRing(QQ,3) sage: G = R.ideal([x^2*y - z, y^2*z - x, z^2*x - y]).groebner_fan() sage: G._gfan_ideal() 'Q[x, y, z]{x^2*y-z,y^2*z-x,x*z^2-y}' """
""" Returns the weight vectors corresponding to the reduced Groebner bases.
EXAMPLES::
sage: r3.<x,y,z> = PolynomialRing(QQ,3) sage: g = r3.ideal([x^3+y,y^3-z,z^2-x]).groebner_fan() sage: g.weight_vectors() [(3, 7, 1), (5, 1, 2), (7, 1, 4), (5, 1, 4), (1, 1, 1), (1, 4, 8), (1, 4, 10)] sage: r4.<x,y,z,w> = PolynomialRing(QQ,4) sage: g4 = r4.ideal([x^3+y,y^3-z,z^2-x,z^3 - w]).groebner_fan() sage: len(g4.weight_vectors()) 23 """
""" Return the multivariate polynomial ring.
EXAMPLES::
sage: R.<x1,x2> = PolynomialRing(QQ,2) sage: gf = R.ideal([x1^3-x2,x2^3-x1-2]).groebner_fan() sage: gf.ring() Multivariate Polynomial Ring in x1, x2 over Rational Field """
""" A string of the reduced Groebner bases of the ideal as output by gfan.
EXAMPLES::
sage: R.<a,b> = PolynomialRing(QQ,2) sage: gf = R.ideal([a^3-b^2,b^2-a-1]).groebner_fan() sage: gf._gfan_reduced_groebner_bases() 'Q[a,b]{{b^6-1+2*b^2-3*b^4,a+1-b^2},{b^2-1-a,a^3-1-a}}' """
""" Return the characteristic of the base ring.
EXAMPLES::
sage: R.<x,y,z> = PolynomialRing(QQ,3) sage: i1 = ideal(x*z + 6*y*z - z^2, x*y + 6*x*z + y*z - z^2, y^2 + x*z + y*z) sage: gf = i1.groebner_fan() sage: gf.characteristic() 0 """
""" EXAMPLES::
sage: R.<x,y,z> = PolynomialRing(QQ, 3, order='lex') sage: G = R.ideal([x^2*y - z, y^2*z - x, z^2*x - y]).groebner_fan() sage: X = G.reduced_groebner_bases() sage: len(X) 33 sage: X[0] [z^15 - z, y - z^11, x - z^9] sage: X[0].ideal() Ideal (z^15 - z, y - z^11, x - z^9) of Multivariate Polynomial Ring in x, y, z over Rational Field sage: X[:5] [[z^15 - z, y - z^11, x - z^9], [-y + z^11, y*z^4 - z, y^2 - z^8, x - z^9], [-y^2 + z^8, y*z^4 - z, y^2*z^3 - y, y^3 - z^5, x - y^2*z], [-y^3 + z^5, y*z^4 - z, y^2*z^3 - y, y^4 - z^2, x - y^2*z], [-y^4 + z^2, y^6*z - y, y^9 - z, x - y^2*z]] sage: R3.<x,y,z> = PolynomialRing(GF(2477),3) sage: gf = R3.ideal([300*x^3-y,y^2-z,z^2-12]).groebner_fan() sage: gf.reduced_groebner_bases() [[z^2 - 12, y^2 - z, x^3 + 933*y], [-y^2 + z, y^4 - 12, x^3 + 933*y], [z^2 - 12, -300*x^3 + y, x^6 - 1062*z], [-828*x^6 + z, -300*x^3 + y, x^12 + 200]] """
""" Return the extra options to the gfan command that are used by this object to account for working modulo a prime or in the presence of extra symmetries.
EXAMPLES::
sage: R.<x,y> = PolynomialRing(QQ,2) sage: gf = R.ideal([x^3-y,y^3-x-1]).groebner_fan() sage: gf._gfan_mod() '' """ #p = self.characteristic() #if p != 0: # mod += ' --mod %s'%p #else: # mod += ''
mod += ' -g'
mod += ' --symmetry'
r""" Returns the gfan output as a string given an input cmd; the default is to produce the list of reduced Groebner bases in gfan format.
EXAMPLES::
sage: R.<x,y> = PolynomialRing(QQ,2) sage: gf = R.ideal([x^3-y,y^3-x-1]).groebner_fan() sage: gf.gfan() 'Q[x,y]\n{{\ny^9-1-y+3*y^3-3*y^6,\nx+1-y^3}\n,\n{\ny^3-1-x,\nx^3-y}\n,\n{\ny-x^3,\nx^9-1-x}\n}\n' """ # todo -- put something in here (?) when self.__symmetry isn't None... raise RuntimeError("Error running gfan command %s on %s"%(cmd, self))
""" Returns an iterator for the reduced Groebner bases.
EXAMPLES::
sage: R.<x,y> = PolynomialRing(QQ,2) sage: gf = R.ideal([x^3-y,y^3-x-1]).groebner_fan() sage: a = gf.__iter__() sage: next(a) [y^9 - 3*y^6 + 3*y^3 - y - 1, -y^3 + x + 1] """
""" Gets a reduced groebner basis
EXAMPLES;
::
sage: R4.<w1,w2,w3,w4> = PolynomialRing(QQ,4) sage: gf = R4.ideal([w1^2-w2,w2^3-1,2*w3-w4^2,w4^2-w1]).groebner_fan() sage: gf[0] [w4^12 - 1, -1/2*w4^2 + w3, -w4^4 + w2, -w4^2 + w1] """
""" Computes and returns a lexicographic reduced Groebner basis for the ideal.
EXAMPLES::
sage: R.<x,y,z> = PolynomialRing(QQ,3) sage: G = R.ideal([x - z^3, y^2 - x + x^2 - z^3*x]).groebner_fan() sage: G.buchberger() [-z^3 + y^2, -z^3 + x] """
""" Returns a polyhedral fan object corresponding to the reduced Groebner bases.
EXAMPLES::
sage: R3.<x,y,z> = PolynomialRing(QQ,3) sage: gf = R3.ideal([x^8-y^4,y^4-z^2,z^2-1]).groebner_fan() sage: pf = gf.polyhedralfan() sage: pf.rays() [[0, 0, 1], [0, 1, 0], [1, 0, 0]] """
""" Return the homogeneity space of a the list of polynomials that define this Groebner fan.
EXAMPLES::
sage: R.<x,y> = PolynomialRing(QQ,2) sage: G = R.ideal([y^3 - x^2, y^2 - 13*x]).groebner_fan() sage: H = G.homogeneity_space() """
""" Render a Groebner fan as sage graphics or save as an xfig file.
More precisely, the output is a drawing of the Groebner fan intersected with a triangle. The corners of the triangle are (1,0,0) to the right, (0,1,0) to the left and (0,0,1) at the top. If there are more than three variables in the ring we extend these coordinates with zeros.
INPUT:
- ``file`` - a filename if you prefer the output saved to a file. This will be in xfig format.
- ``shift`` - shift the positions of the variables in the drawing. For example, with shift=1, the corners will be b (right), c (left), and d (top). The shifting is done modulo the number of variables in the polynomial ring. The default is 0.
- ``larger`` - bool (default: ``False``); if ``True``, make the triangle larger so that the shape of the Groebner region appears. Affects the xfig file but probably not the sage graphics (?)
- ``rgbcolor`` - This will not affect the saved xfig file, only the sage graphics produced.
- ``polyfill`` - Whether or not to fill the cones with a color determined by the highest degree in each reduced Groebner basis for that cone.
- ``scale_colors`` - if True, this will normalize color values to try to maximize the range
EXAMPLES::
sage: R.<x,y,z> = PolynomialRing(QQ,3) sage: G = R.ideal([y^3 - x^2, y^2 - 13*x,z]).groebner_fan() sage: test_render = G.render()
::
sage: R.<x,y,z> = PolynomialRing(QQ,3) sage: G = R.ideal([x^2*y - z, y^2*z - x, z^2*x - y]).groebner_fan() sage: test_render = G.render(larger=True)
TESTS:
Testing the case where the number of generators is < 3. Currently, this should raise a ``NotImplementedError`` error.
::
sage: R.<x,y> = PolynomialRing(QQ, 2) sage: R.ideal([y^3 - x^2, y^2 - 13*x]).groebner_fan().render() Traceback (most recent call last): ... NotImplementedError """ cmd += ' --shiftVariables %s'%shift open(file,'w').write(s) if scale_colors: vmins = [min([q[i] for q in vals]) for i in (0,1,2)] vmaxs = [max([q[i] for q in vals]) for i in (0,1,2)] for i in (0,1,2): if vmaxs[i] == vmins[i]: vmaxs[i] = vmins[i] + .01 for index in range(len(sp3)): col = [1-(vals[index][i]-vmins[i])/(vmaxs[i]-vmins[i]) for i in (0,1,2)] r_lines = r_lines + polygon(sp3[index], rgbcolor = col) else: for index in range(len(sp3)): r_lines = r_lines + polygon(sp3[index], rgbcolor = vals[i]) else: vmax = vmin + .01 else: for index in range(len(sp3)): r_lines = r_lines + polygon(sp3[index], hue = vals[i])
""" A simple utility function for converting a facet normal to an inequality form.
EXAMPLES::
sage: R.<x,y> = PolynomialRing(QQ,2) # dummy stuff to get a gfan object sage: gf = R.ideal([x^2+y,y^2]).groebner_fan() sage: gf._cone_to_ieq([[1,2,3,4]]) [[0, 1, 2, 3, 4]] """
""" Takes a 4-d vector and projects it onto the plane perpendicular to (1,1,1,1). Stretches by a factor of 2 as well, since this is only for graphical display purposes.
INPUT:
- ``fpoint`` - a list of four numbers
EXAMPLES::
sage: R.<x> = PolynomialRing(QQ,1) # dummy stuff to get a gfan object sage: gf = R.ideal([x^2]).groebner_fan() sage: gf._embed_tetra([1/2,1/2,1/2,1/2]) [0, 0, 0] """
""" A utility function that takes a list of 4d polytopes, projects them to 3d, and returns a list of edges.
INPUT:
- ``polyhedral_data`` -- an object with 4d vertex and adjacency information
OUTPUT:
- ``edges`` -- a list of edges in 3d - each list item is a pair of points
EXAMPLES::
sage: R4.<w,x,y,z> = PolynomialRing(QQ,4) sage: gf = R4.ideal([w^2-x,x^2-y,y^2-z,z^2-1]).groebner_fan() sage: g_cone = gf[0].groebner_cone() sage: g_cone_facets = g_cone.facets() sage: g_cone_ieqs = gf._cone_to_ieq(g_cone_facets) sage: cone_data = Polyhedron(ieqs = g_cone_ieqs, eqns = [[1,-1,-1,-1,-1]]) sage: cone_lines = gf._4d_to_3d(cone_data) sage: cone_lines [[[1, 1, -1], [1, -1/3, 1/3]], [[1, 1, -1], [-1/7, 3/7, 5/7]], [[1, 1, -1], [-3/5, -1/3, -1/5]], [[1, -1/3, 1/3], [-1/7, 3/7, 5/7]], [[1, -1/3, 1/3], [-3/5, -1/3, -1/5]], [[-1/7, 3/7, 5/7], [-3/5, -1/3, -1/5]]] """ except Exception: print(adj) print('tpoints: ' + str(tpoints)) print('fpoints: ' + str(fpoints)) print(adjacent_vertex) print(polyhedral_data.ieqs()) raise RuntimeError(adj)
""" For a Groebner fan of an ideal in a ring with four variables, this function intersects the fan with the standard simplex perpendicular to (1,1,1,1), creating a 3d polytope, which is then projected into 3 dimensions. The edges of this projected polytope are returned as lines.
EXAMPLES::
sage: R4.<w,x,y,z> = PolynomialRing(QQ,4) sage: gf = R4.ideal([w^2-x,x^2-y,y^2-z,z^2-x]).groebner_fan() sage: three_d = gf.render3d()
TESTS:
Now test the case where the number of generators is not 4. Currently, this should raise a ``NotImplementedError`` error.
::
sage: P.<a,b,c> = PolynomialRing(QQ, 3, order="lex") sage: sage.rings.ideal.Katsura(P, 3).groebner_fan().render3d() Traceback (most recent call last): ... NotImplementedError """ # Now the cones are intersected with a plane: #This is really just for debugging for x in cone_info: print(x.inequalities() + ([1,1,0,0,0],[1,0,1,0,0],[1,0,0,1,0],[1,0,0,0,1])) print(x.equations()) print() ([1,1,0,0,0],[1,0,1,0,0],[1,0,0,1,0],[1,0,0,0,1]), eqns = x.equations()) for x in cone_info] except Exception: print(cone_data._rays) raise RuntimeError
""" Return various statistics about this Groebner fan.
EXAMPLES::
sage: R.<x,y> = PolynomialRing(QQ) sage: G = R.ideal([y^3 - x^2, y^2 - 13*x]).groebner_fan() sage: G._gfan_stats() {'Dimension of homogeneity space': 0, 'Maximal number of polynomials in Groebner basis': 3, 'Maximal number of terms in Groebner basis': 6, 'Maximal total degree of a Groebner basis': 4, 'Minimal total degree of a Groebner basis': 2, 'Number of reduced Groebner bases': 3, 'Number of variables': 2} """
""" Return the dimension of the homogeneity space.
EXAMPLES::
sage: R.<x,y> = PolynomialRing(QQ,2) sage: G = R.ideal([y^3 - x^2, y^2 - 13*x]).groebner_fan() sage: G.dimension_of_homogeneity_space() 0 """
""" Return the maximal total degree of any Groebner basis.
EXAMPLES::
sage: R.<x,y> = PolynomialRing(QQ,2) sage: G = R.ideal([y^3 - x^2, y^2 - 13*x]).groebner_fan() sage: G.maximal_total_degree_of_a_groebner_basis() 4 """
""" Return the minimal total degree of any Groebner basis.
EXAMPLES::
sage: R.<x,y> = PolynomialRing(QQ,2) sage: G = R.ideal([y^3 - x^2, y^2 - 13*x]).groebner_fan() sage: G.minimal_total_degree_of_a_groebner_basis() 2 """
""" Return the number of reduced Groebner bases.
EXAMPLES::
sage: R.<x,y> = PolynomialRing(QQ,2) sage: G = R.ideal([y^3 - x^2, y^2 - 13*x]).groebner_fan() sage: G.number_of_reduced_groebner_bases() 3 """
""" Return the number of variables.
EXAMPLES::
sage: R.<x,y> = PolynomialRing(QQ,2) sage: G = R.ideal([y^3 - x^2, y^2 - 13*x]).groebner_fan() sage: G.number_of_variables() 2
::
sage: R = PolynomialRing(QQ,'x',10) sage: R.inject_variables(globals()) Defining x0, x1, x2, x3, x4, x5, x6, x7, x8, x9 sage: G = ideal([x0 - x9, sum(R.gens())]).groebner_fan() sage: G.number_of_variables() 10 """
""" Return a tropical basis for the tropical curve associated to this ideal.
INPUT:
- ``check`` - bool (default: True); if True raises a ValueError exception if this ideal does not define a tropical curve (i.e., the condition that R/I has dimension equal to 1 + the dimension of the homogeneity space is not satisfied).
EXAMPLES::
sage: R.<x,y,z> = PolynomialRing(QQ,3, order='lex') sage: G = R.ideal([y^3-3*x^2, z^3-x-y-2*y^3+2*x^2]).groebner_fan() sage: G Groebner fan of the ideal: Ideal (-3*x^2 + y^3, 2*x^2 - x - 2*y^3 - y + z^3) of Multivariate Polynomial Ring in x, y, z over Rational Field sage: G.tropical_basis() [-3*x^2 + y^3, 2*x^2 - x - 2*y^3 - y + z^3, 3/4*x + y^3 + 3/4*y - 3/4*z^3] """
raise ValueError("The ideal does not define a tropical curve.")
""" See the documentation for self[0].interactive(). This does not work with the notebook.
EXAMPLES::
sage: print("This is not easily doc-testable; please write a good one!") This is not easily doc-testable; please write a good one! """ self[0].interactive(*args, **kwds)
""" Returns information about the tropical intersection of the polynomials defining the ideal. This is the common refinement of the outward-pointing normal fans of the Newton polytopes of the generators of the ideal. Note that some people use the inward-pointing normal fans.
INPUT:
- ``parameters`` (optional) - a list of variables to be considered as parameters - ``symmetry_generators`` (optional) - generators of the symmetry group
OUTPUT: a TropicalPrevariety object
EXAMPLES::
sage: R.<x,y,z> = PolynomialRing(QQ,3) sage: I = R.ideal(x*z + 6*y*z - z^2, x*y + 6*x*z + y*z - z^2, y^2 + x*z + y*z) sage: gf = I.groebner_fan() sage: pf = gf.tropical_intersection() sage: pf.rays() [[-2, 1, 1]]
sage: R.<x,y,z> = PolynomialRing(QQ,3) sage: f1 = x*y*z - 1 sage: f2 = f1*(x^2 + y^2 + z^2) sage: f3 = f2*(x + y + z - 1) sage: I = R.ideal([f1,f2,f3]) sage: gf = I.groebner_fan() sage: pf = gf.tropical_intersection(symmetry_generators = '(1,2,0),(1,0,2)') sage: pf.rays() [[-2, 1, 1], [1, -2, 1], [1, 1, -2]]
sage: R.<x,y,z> = QQ[] sage: I = R.ideal([(x+y+z)^2-1,(x+y+z)-x,(x+y+z)-3]) sage: GF = I.groebner_fan() sage: TI = GF.tropical_intersection() sage: TI.rays() [[-1, 0, 0], [0, -1, -1], [1, 1, 1]] sage: GF = I.groebner_fan() sage: TI = GF.tropical_intersection(parameters=[y]) sage: TI.rays() [[-1, 0, 0]] """
""" Returns the mixed volume of the generators of this ideal (i.e. this is not really an ideal property, it can depend on the generators used). The generators must give a square system (as many polynomials as variables).
EXAMPLES::
sage: R.<x,y,z> = QQ[] sage: example_ideal = R.ideal([x^2-y-1,y^2-z-1,z^2-x-1]) sage: gf = example_ideal.groebner_fan() sage: mv = gf.mixed_volume() sage: mv 8
sage: R2.<x,y> = QQ[] sage: g1 = 1 - x + x^7*y^3 + 2*x^8*y^4 sage: g2 = 2 + y + 3*x^7*y^3 + x^8*y^4 sage: example2 = R2.ideal([g1,g2]) sage: example2.groebner_fan().mixed_volume() 15 """ raise UserWarning('not a square system')
""" A class for representing reduced Groebner bases as produced by gfan.
INPUT:
- ``groebner_fan`` - a GroebnerFan object from an ideal
- ``gens`` - the generators of the ideal
- ``gfan_gens`` - the generators as a gfan string
EXAMPLES::
sage: R.<a,b> = PolynomialRing(QQ,2) sage: gf = R.ideal([a^2-b^2,b-a-1]).groebner_fan() sage: from sage.rings.polynomial.groebner_fan import ReducedGroebnerBasis sage: ReducedGroebnerBasis(gf,gf[0],gf[0]._gfan_gens()) [b - 1/2, a + 1/2] """
""" Returns the reduced Groebner basis as a string.
EXAMPLES::
sage: R.<z1,zz1> = PolynomialRing(QQ,2) sage: gf = R.ideal([z1^2*zz1-1,zz1-2]).groebner_fan() sage: rgb1 = gf.reduced_groebner_bases()[0] sage: rgb1 # indirect doctest [zz1 - 2, z1^2 - 1/2] """
""" Returns the reduced Groebner basis as a string in gfan format.
EXAMPLES::
sage: R.<z1,zz1> = PolynomialRing(QQ,2) sage: gf = R.ideal([z1^2*zz1-1,zz1-2]).groebner_fan() sage: rgb1 = gf.reduced_groebner_bases()[0] sage: rgb1._gfan_gens() '{zz1-2,z1^2-1/2}' """
""" Returns a description of the Groebner fan this basis was derived from.
EXAMPLES::
sage: R.<z1,zz1> = PolynomialRing(QQ,2) sage: gf = R.ideal([z1^2*zz1-1,zz1-2]).groebner_fan() sage: rgb1 = gf.reduced_groebner_bases()[0] sage: rgb1._gfan() Groebner fan of the ideal: Ideal (z1^2*zz1 - 1, zz1 - 2) of Multivariate Polynomial Ring in z1, zz1 over Rational Field """
inequalities=False, weight=False): """ Do an interactive walk of the Groebner fan starting at this reduced Groebner basis.
EXAMPLES::
sage: R.<x,y> = PolynomialRing(QQ,2) sage: G = R.ideal([y^3 - x^2, y^2 - 13*x]).groebner_fan() sage: G[0].interactive() # not tested Initializing gfan interactive mode ********************************************* * Press control-C to return to Sage * ********************************************* .... """ cmd = 'gfan_interactive' if latex: cmd += ' -L' if flippable: cmd += ' -f' if wall: cmd += ' -w' if inequalities: cmd += ' -i' if weight: cmd += ' -W' cmd += self.__groebner_fan._gfan_mod() E = pexpect.spawn(cmd) print("Initializing gfan interactive mode") #E.sendline(self._gfan_ideal()) E.sendline(self.__gfan_gens) print("*" * 45) print("* Press control-C to return to Sage *") print("*" * 45) try: E.interact() except OSError: print("Returning to Sage.")
""" Return defining inequalities for the full-dimensional Groebner cone associated to this marked minimal reduced Groebner basis.
INPUT:
- ``restrict`` - bool (default: False); if True, add an inequality for each coordinate, so that the cone is restricted to the positive orthant.
OUTPUT: tuple of integer vectors
EXAMPLES::
sage: R.<x,y> = PolynomialRing(QQ,2) sage: G = R.ideal([y^3 - x^2, y^2 - 13*x]).groebner_fan() sage: poly_cone = G[1].groebner_cone() sage: poly_cone.facets() [[-1, 2], [1, -1]] sage: [g.groebner_cone().facets() for g in G] [[[0, 1], [1, -2]], [[-1, 2], [1, -1]], [[-1, 1], [1, 0]]] sage: G[1].groebner_cone(restrict=True).facets() [[-1, 2], [1, -1]] """
""" Return the ideal generated by this basis.
EXAMPLES::
sage: R.<x,y,z> = PolynomialRing(QQ,3) sage: G = R.ideal([x - z^3, y^2 - 13*x]).groebner_fan() sage: G[0].ideal() Ideal (-13*z^3 + y^2, -z^3 + x) of Multivariate Polynomial Ring in x, y, z over Rational Field """
|