Ticket #4539: letterplace.py

File letterplace.py, 3.8 KB (added by PolyBoRi, 5 years ago)
Line 
1from sage.libs.singular.letterplace import temporary_MPolynomialRing_from_letterplace_ring
2from sage.libs.singular.function import singular_function
3from sage.algebras.free_algebra import FreeAlgebra
4from sage.libs.singular.function import lib, singular_function
5from sage.rings.polynomial.polynomial_ring_constructor import PolynomialRing
6
7def iter_free_monomial(monomial, gens):
8    for (index, exponent) in monomial:
9        for i in xrange(exponent):
10            yield gens[index]
11
12def freegb(ideal, deg_bound):
13    """
14    sage: F.<x,y,z> = FreeAlgebra(QQ, 3); F
15    Free Algebra on 3 generators (x, y, z) over Rational Field
16    sage: l=[2*x*z*x+y*x*y, 3*x*y+x*z]
17    sage: freegb(l, 10)
18    [3*y*x*z^7*y + y*x*z^8, 3*y*x*z^6*y + y*x*z^7, y*x*z^6*x*z + 314928*y^2*x*z^2*x^5, 3*y*x*z^5*y + y*x*z^6, y*x*z^5*x*z - 17496*y^2*x*z^2*x^4, 3*y*x*z^4*y + y*x*z^5, y*x*z^4*x*z + 972*y^2*x*z^2*x^3, 3*y*x*z^4*x^2*z*y + y*x*z^4*x^2*z^2, 3*y*x*z^3*y + y*x*z^4, y*x*z^3*x*z - 54*y^2*x*z^2*x^2, 3*y*x*z^3*x^2*z^2*y + y*x*z^3*x^2*z^3, 3*y*x*z^3*x^2*z*y + y*x*z^3*x^2*z^2, 3*y*x*z^3*x^3*z*y + y*x*z^3*x^3*z^2, 3*y*x*z^2*y + y*x*z^3, y*x*z^2*x*z + 3*y^2*x*z^2*x, 3*y*x*z^2*x^2*z^3*y + y*x*z^2*x^2*z^4, 3*y*x*z^2*x^2*z^2*y + y*x*z^2*x^2*z^3, y*x*z^2*x^2*z^2*x*z + 3*y^2*x*z^2*x^2*z^2*x, 3*y*x*z^2*x^2*z*y + y*x*z^2*x^2*z^2, 3*y*x*z^2*x^3*z^2*y + y*x*z^2*x^3*z^3, 3*y*x*z^2*x^3*z*y + y*x*z^2*x^3*z^2, 3*y*x*z^2*x^4*z*y + y*x*z^2*x^4*z^2, 3*y*x*z*y + y*x*z^2, x*z*y^6*x*z - 7776*y*x*z^2*x^6, x*z*y^5*x*z - 1296*y*x*z^2*x^5, x*z*y^4*x*z - 216*y*x*z^2*x^4, x*z*y^3*x*z - 36*y*x*z^2*x^3, x*z*y^2*x*z - 6*y*x*z^2*x^2, x*z*y*x*z - y*x*z^2*x, 6*x*z*x - y*x*z, 3*x*y + x*z]
19    """
20    lib("freegb.lib")
21    free_algebra=ideal[0].parent()
22    base_ring = free_algebra.base_ring()
23    gens = free_algebra.gens()
24    variable_names = [str(g) for g in gens]
25    pre_letter_place_ring=PolynomialRing(base_ring, variable_names)
26    make_letter_place_ring = singular_function("makeLetterplaceRing")
27    ring_wrap = make_letter_place_ring(10, ring=pre_letter_place_ring)
28   
29    (to_letter_place, from_letter_place) = \
30        temporary_MPolynomialRing_from_letterplace_ring(ring_wrap, 
31        base_ring, gens,  deg_bound
32        )
33   
34    singular_ring = iter(from_letter_place.keys()).next().parent()
35   
36    polynomial_list=[]
37    for p in ideal:
38        sum_p = 0
39        for term in p:
40            (coef, monomial) = term
41            m = 1
42            for (i,v) in enumerate(iter_free_monomial(monomial, gens)):
43                m=m*to_letter_place[v][i]
44            m=m*coef
45            sum_p = sum_p+m
46        polynomial_list.append(sum_p)
47   
48   
49    free_gbasis=singular_function("freeGBasis")
50    from sage.libs.singular.option import LibSingularOptions
51    libsingular_options = LibSingularOptions()
52    bck = int(libsingular_options)
53    #letter place needs these options
54    libsingular_options['redTail'] = True
55    libsingular_options['redSB'] = True
56    libsingular_options(bck)
57   
58    singular_system=singular_function("system")
59    gb = singular_system("freegb",
60        singular_ring.ideal(polynomial_list), 
61        deg_bound, 
62        len(gens))
63
64    result = []
65    def sort_key(index_generator_tuple):
66        return index_generator_tuple[0]
67    for p in gb:
68        sum_p=free_algebra(0)
69        for (coef, monomial) in p:
70            m=free_algebra(1)
71            index_generator_list=[]
72            for v in monomial.variables():
73                #there exist no exponent > 1 in the result freegb
74                (index, alg_gen) = from_letter_place[v]
75                index_generator_list.append((index, alg_gen))
76            index_generator_list = sorted(index_generator_list, key=sort_key)
77           
78            for (index, alg_gen) in index_generator_list:
79                m = m * alg_gen
80            m=coef * m
81            sum_p = sum_p + m
82        result.append(sum_p)
83    return result