API tools faq
paste
Guest User

Untitled

a guest
Feb 21st, 2018
56
0
Never
Add comment
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
text 48.29 KB | None | 0 0
raw download clone embed print report
  1. r"""
  2.  
  3. Educational Versions of the F5 Algorithm to compute tropical Groebner bases.
  4.  
  5. We have followed [F02], [AP11] [EF17] and [VY17].
  6.  
  7. No attempt was made to optimize the algorithm as the emphasis of
  8. these implementations is to obtain an as clean and easy presentation as possible. To compute a
  9. Groebner basis in Sage efficiently use the
  10. :meth:`sage.rings.polynomial.multi_polynomial_ideal.MPolynomialIdeal.groebner_basis()`
  11. method on multivariate polynomial objects.
  12.  
  13. The main algorithm is ``tentative_F5Trop``; secondary algorithms of importance are ``SymbolicPreprocessing``,
  14. ``RowEchelonMac_and_listpivots`` and ``UpdateGpaires_trop``.
  15.  
  16. The idea of the algorithm is to use an adapted Buchberger termination criterion
  17. on $S$-pairs while reducing all computations to linear algebra
  18. and eliminating useless polynomials with the F5 elimination criterion.
  19. The core part of the main algorithm is a while loop
  20. while there are still S-polynomials, where it applies:
  21.  
  22. * ``SymbolicPreprocessing`` to construct a matrix with the polynomials of the pairs of a given degree d and reducers to reduce them. It uses the F5 elimination criterion.
  23. * ``RowEchelonMac_and_listpivots`` to echelonize the previous matrix.
  24. * ``UpdateGpaires_trop`` to update the list in construction of a Groebner basis (more precisely, an $\mathfrak{S}$-Groebner basis), and obtain the new S-polynomials.
  25.  
  26.  
  27. EXAMPLES:
  28.  
  29. Consider Katsura-5 over QQ w.r.t. 2-adic valuation and a ``degrevlex`` ordering.::
  30.  
  31. S.<a,b,c,d,e>=QQ[]
  32. P.<x,y,z,s,t>=Qp(2,50)[]
  33. w=[1,-2,4,-8,16]
  34. listpol = sage.rings.ideal.Katsura(S).gens()
  35.  
  36. G5=reduced_trop_GB(listpol,P,w)
  37.  
  38. REFERENCES:
  39.  
  40. .. [F02] Jean-Charles Faugère. *A new efficient algorithm for computing Gröbner bases without reduction to zero (F5)*. In Proceedings of the 2002 international symposium on Symbolic and algebraic computation, ISSAC '02, pages 75-83, New York, NY, USA, 2002. ACM.
  41. .. [AP11] Alberto Arri and John Perry. *The F5 criterion revised*. In Journal of Symbolic Computation, 2011. Corrigendum in 2017.
  42. .. [EF17] Christian Eder and Jean-Charles Faugère. *A survey on signature-based algorithms for computing Gröbner bases*. In Journal of Symbolic Computation, 2017.
  43. .. [VY17] T. Vaccon, K.Yokoyama, A Tropical F5 Algorithm, proceedings of ISSAC 2017.
  44.  
  45. AUTHOR:
  46.  
  47. - Tristan Vaccon (2016--2018)
  48. """
  49.  
  50. from copy import copy
  51. from sage.structure.sequence import Sequence
  52. from sage.matrix.matrix_space import MatrixSpace
  53. from sage.matrix.constructor import matrix, Matrix
  54. from sage.rings.polynomial.multi_polynomial_element import MPolynomial_polydict
  55. from sage.rings.polynomial.polydict import PolyDict, ETuple
  56.  
  57. def Make_Macaulay_Matrix_sparse_dict(list_poly, monom):
  58. r"""
  59.  
  60. Build the Macaulay matrix for the polynomials of degree d in list_poly. It is represented as a couple of a matrix and the list of the signatures of its rows, assuming compatibility.
  61. The matrix is sparse.
  62.  
  63. INPUT:
  64.  
  65. - ``list_poly`` -- a list of degree d polynomials
  66. - ``monom`` -- the list of the monomials of degree d, ordered decreasingly
  67.  
  68. OUTPUT:
  69.  
  70. The Macaulay matrix defined by list_poly
  71.  
  72. EXAMPLES::
  73.  
  74. sage: S.<x,y,z>=QQ[]
  75. sage: from sage.rings.polynomial.padics.F5Trop_doctest import Make_Macaulay_Matrix_sparse_dict
  76. sage: list_poly = [[0,1,x**2+y**2], [1,1,x*y], [2,1,y*z], [3,1,z**2]]
  77. sage: Mac = Make_Macaulay_Matrix_sparse_dict(list_poly,[x**2, x*y, y**2, x*z, y*z, z**2])
  78. sage: Mac[0]
  79. [1 0 1 0 0 0]
  80. [0 1 0 0 0 0]
  81. [0 0 0 0 1 0]
  82. [0 0 0 0 0 1]
  83. sage: Mac[1]
  84. [[0, 1], [1, 1], [2, 1], [3, 1]]
  85.  
  86. """
  87.  
  88. R = list_poly[0][2].parent()
  89. d = list_poly[0][2].degree()
  90. l = len(monom)
  91.  
  92. dict_monom = {monom[aa]:aa for aa in range(l)}
  93.  
  94.  
  95. nrows = len(list_poly)
  96.  
  97. listsign = []
  98. Mac = MatrixSpace(R.base_ring(),nrows,l,sparse=True)(0)
  99. for u in range(nrows):
  100. fuple=list_poly[u]
  101. listsign.append([fuple[0],fuple[1]])
  102. f = fuple[2]
  103. list = f.monomials()
  104. for mon in list:
  105. if dict_monom.has_key(mon):
  106. j = dict_monom[mon]
  107. Mac[u,j] = f.monomial_coefficient(mon)
  108. Mac = [Mac, listsign]
  109. return Mac
  110.  
  111.  
  112. def Reconstruction(Mactilde,i,monom,sign):
  113. r"""
  114.  
  115. Computes the polynomial with signature, corresponding to the row of Mactilde of index i
  116.  
  117.  
  118. INPUT:
  119.  
  120. - ``Mactilde`` -- a Macaulay matrix (just the matrix)
  121. - ``i`` -- an integer
  122. - ``monom`` -- the list of the monomials corresponding to the columns of Mactilde
  123. - ``sign`` -- the signature corresponding to the row of Mactilde of index i
  124.  
  125. OUTPUT:
  126.  
  127. a polynomial with signature, corresponding to the row of Mactilde of index i
  128.  
  129. EXAMPLES::
  130.  
  131. sage: from sage.rings.polynomial.padics.F5Trop_doctest import*
  132. sage: S.<x,y,z>=QQ[]
  133. sage: mat = Matrix(QQ,2,3,[0,1,2,1,-1,-2])
  134. sage: monom = [x,y,z]
  135. sage: sign = [0,1]
  136. sage: Reconstruction(mat,1,monom,sign)
  137. [0, 1, x - y - 2*z]
  138. """
  139. r = Mactilde.row(i)
  140. m = Mactilde.ncols()
  141. pol = sum(r[j]*monom[j] for j in range(m))
  142. return sign+[pol]
  143.  
  144.  
  145.  
  146. def Make_Macaulay_Matrix_sparse_dict_no_sign(list_poly, monom):
  147. r"""
  148.  
  149. Make the sparse Macaulay matrix for the polynomials of degree d in list_poly
  150.  
  151. INPUT:
  152.  
  153. - ``list_poly`` -- a list of degree d polynomials
  154. - ``monom`` -- the list of the monomials of degree d, ordered decreasingly
  155.  
  156. OUTPUT:
  157.  
  158. The Macaulay matrix defined by list_poly
  159.  
  160. EXAMPLES::
  161.  
  162. sage: S.<x,y,z>=QQ[]
  163. sage: from sage.rings.polynomial.padics.F5Trop_doctest import Make_Macaulay_Matrix_sparse_dict_no_sign
  164. sage: list_poly = [x^2+y^2,x*y,y*z,z^2]
  165. sage: Mac = Make_Macaulay_Matrix_sparse_dict_no_sign(list_poly,[x^2, x*y, y^2, x*z, y*z, z^2])
  166. sage: Mac
  167. [1 0 1 0 0 0]
  168. [0 1 0 0 0 0]
  169. [0 0 0 0 1 0]
  170. [0 0 0 0 0 1]
  171.  
  172.  
  173. """
  174.  
  175. R = list_poly[0].parent()
  176. d = list_poly[0].degree()
  177. l = len(monom)
  178.  
  179. dict_monom = {monom[aa]:aa for aa in range(l)}
  180.  
  181.  
  182. nrows = len(list_poly)
  183.  
  184.  
  185. Mac = MatrixSpace(R.base_ring(),nrows,l,sparse=True)(0)
  186. for u in range(nrows):
  187. f = list_poly[u]
  188. list = f.monomials()
  189. for mon in list:
  190. if dict_monom.has_key(mon):
  191. j = dict_monom[mon]
  192. Mac[u,j] = f.monomial_coefficient(mon)
  193. return Mac
  194.  
  195. def LT_trop(f,R,w):
  196. r"""
  197.  
  198. Computes the leading term of f according to the degree-refining tropical term order
  199. given by the monomial ordering on R and the weight w
  200.  
  201. INPUT:
  202.  
  203. - ``f`` -- a polynomial
  204. - ``R`` -- a polynomial ring of dimension
  205. - ``w`` -- a list of n rational numbers
  206.  
  207. OUTPUT:
  208.  
  209. the leading term of f
  210.  
  211. EXAMPLES::
  212.  
  213. sage: S.<x,y,z>=QQ[]
  214. sage: from sage.rings.polynomial.padics.F5Trop_doctest import LT_trop
  215. sage: f = 2*x^2+3*y^2+z^2+1
  216. sage: R1.<x1,y1,z1> = Qp(2)[]
  217. sage: R2 = PolynomialRing(Qp(3),['x2','y2','z2'] , order = 'lex')
  218. sage: w1 = [0,0,0]
  219. sage: w2 = [0,0,-1]
  220. sage: LT_trop(f,R1,w1)
  221. 3*y^2
  222. sage: LT_trop(f,R2,w2)
  223. z^2
  224. """
  225. d=f.degree()
  226. lmon = [u for u in f.monomials() if u.degree()==d]
  227. lcoeff =[f.monomial_coefficient(u) for u in lmon]
  228. nvar = len(f.parent().gens())
  229. K = R.base_ring()
  230. def to_R(mon):
  231. d = PolyDict({ETuple(mon.exponents()[0]): K.one()}, force_etuples=False)
  232. return MPolynomial_polydict(R, d)
  233.  
  234. if len(lmon) == 0:
  235. return f
  236. else:
  237. imin = 0
  238. mon = lmon[0]
  239. trop_weight0 = K(lcoeff[0]).valuation()+sum([lmon[0].degrees()[vv]*w[vv] for vv in range(nvar)])
  240. for uu in range(len(lmon)):
  241. trop_weight = K(lcoeff[uu]).valuation()+sum([lmon[uu].degrees()[ww]*w[ww] for ww in range(nvar)])
  242. if trop_weight < trop_weight0 or (trop_weight == trop_weight0 and to_R(lmon[uu])> to_R(mon)):
  243. imin = uu
  244. mon = lmon[uu]
  245. trop_weight0 = trop_weight
  246. return lcoeff[imin]*lmon[imin]
  247.  
  248. def LM_trop(f,R,w):
  249. r"""
  250.  
  251. Computes the monomial of the leading term of f according to the degree-refining tropical term order
  252. given by the monomial ordering on R and the weight w
  253.  
  254. INPUT:
  255.  
  256. - ``f`` -- a polynomial
  257. - ``R`` -- a polynomial ring of dimension
  258. - ``w`` -- a list of n rational numbers
  259.  
  260. OUTPUT:
  261.  
  262. a monomial of f, the monomial of the leading term of f
  263.  
  264. EXAMPLES::
  265.  
  266. sage: S.<x,y,z>=QQ[]
  267. sage: from sage.rings.polynomial.padics.F5Trop_doctest import LM_trop
  268. sage: f = 2*x^2+3*y^2+z^2+1
  269. sage: R1.<x1,y1,z1> = Qp(2)[]
  270. sage: R2 = PolynomialRing(Qp(3),['x2','y2','z2'] , order = 'lex')
  271. sage: w1 = [0,0,0]
  272. sage: w2 = [0,0,-1]
  273. sage: LM_trop(f,R1,w1)
  274. y^2
  275. sage: LM_trop(f,R2,w2)
  276. z^2
  277. """
  278. d=f.degree()
  279. lmon = [u for u in f.monomials() if u.degree()==d]
  280. lcoeff =[f.monomial_coefficient(u) for u in lmon]
  281. nvar = len(f.parent().gens())
  282. K = R.base_ring()
  283. def to_R(mon):
  284. d = PolyDict({ETuple(mon.exponents()[0]): K.one()}, force_etuples=False)
  285. return MPolynomial_polydict(R, d)
  286.  
  287. if len(lmon) == 0:
  288. return f
  289. else:
  290. imin = 0
  291. mon = lmon[0]
  292. trop_weight0 = K(lcoeff[0]).valuation()+sum([lmon[0].degrees()[vv]*w[vv] for vv in range(nvar)])
  293. for uu in range(len(lmon)):
  294. trop_weight = K(lcoeff[uu]).valuation()+sum([lmon[uu].degrees()[ww]*w[ww] for ww in range(nvar)])
  295. if trop_weight < trop_weight0 or (trop_weight == trop_weight0 and to_R(lmon[uu])> to_R(mon)):
  296. imin = uu
  297. mon = lmon[uu]
  298. trop_weight0 = trop_weight
  299. return lmon[imin]
  300.  
  301.  
  302. def F5_criterion_trop_affine(sign,G,R,w, list_degrees):
  303. r"""
  304.  
  305. Inside the F5 algorithm
  306. Check the F5 criterion for the signature sign = (i, x^alpha) and G
  307. an S-GB in construction.
  308.  
  309. INPUT:
  310.  
  311. - ``sign`` -- a (potential) signature of a polynomial
  312. - ``G`` -- a list of polynomials with signature
  313. - ``R`` -- a polynomial ring
  314. - ``w`` -- a list of n rational numbers
  315. - ``list_degrees`` -- a list of integers, representing the degrees of the polynomials in F
  316.  
  317. OUTPUT:
  318.  
  319. a boolean, testing whether the F5 criterion for sign and G is passed
  320.  
  321. EXAMPLES::
  322.  
  323. sage: S.<x,y,z>=QQ[]
  324. sage: from sage.rings.polynomial.padics.F5Trop_doctest import *
  325. sage: sign = [1,x^3]
  326. sage: G = [[0,S(1),x^2], [1,S(1),x*y], [2,S(1),y*z], [3,S(1),2*x^2+3*y^2+z^2]]
  327. sage: R1.<x1,y1,z1> = Qp(2)[]
  328. sage: R2 = PolynomialRing(Qp(3),['x2','y2','z2'] , order = 'lex')
  329. sage: w1 = [0,0,0]
  330. sage: w2 = [0,0,-1]
  331. sage: F5_criterion_trop_affine(sign,G, R1,w1,[2,2,2,2])
  332. True
  333. sage: sign2 = [0,x^2]
  334. sage: sign3 = [4,y^2]
  335. sage: F5_criterion_trop_affine(sign2,G, R2, w2,[2,2,2,2])
  336. False
  337. sage: F5_criterion_trop_affine(sign3,G, R1, w1,[2,2,2,2,4])
  338. True
  339. sage: F5_criterion_trop_affine(sign3,G, R2, w2,[2,2,2,2,4])
  340. False
  341.  
  342. """
  343. P = G[0][2].parent()
  344. for g in G:
  345. ei = g[0]
  346. glm = LM_trop(g[2],R,w)
  347. if ei < sign[0]:
  348. if P.monomial_divides(glm,sign[1]):
  349. quot = P.monomial_quotient(sign[1],glm)
  350. if sign[1].degree() >= quot.degree()+g[1].degree()+list_degrees[g[0]]:
  351. return True
  352. return False
  353.  
  354. def Reducible_large_resp_sign_trop(f,g,R,w):
  355. r"""
  356.  
  357. Inside the F5 algorithm, f and g are of the form [i,x^alpha, poly].
  358. Tests whether g can reduce f while respecting their signatures with non-strict inequality
  359.  
  360. INPUT:
  361.  
  362. - ``f`` -- a polynomial with signature
  363. - ``g`` -- a polynomial with signature
  364. - ``R`` -- a polynomial ring of dimension
  365. - ``w`` -- a list of n rational numbers
  366.  
  367. OUTPUT:
  368.  
  369. a boolean, testing whether g can reduce f while respecting their signatures
  370.  
  371. EXAMPLES::
  372.  
  373. sage: S.<x,y,z>=QQ[]
  374. sage: R1.<x1,y1,z1> = Qp(2)[]
  375. sage: R2 = PolynomialRing(Qp(3),['x2','y2','z2'] , order = 'lex')
  376. sage: w1 = [0,0,0]
  377. sage: w2 = [0,0,-1]
  378. sage: from sage.rings.polynomial.padics.F5Trop_doctest import Reducible_large_resp_sign_trop
  379. sage: g = [0,x,x^2]
  380. sage: f = [1,x^2,x^3*y]
  381. sage: f2 = [0,x^2,x^3]
  382. sage: f3 = [0,S(1),x^3]
  383. sage: Reducible_large_resp_sign_trop(f,g,R1,w1)
  384. False
  385. sage: Reducible_large_resp_sign_trop(f2,g,R1,w1)
  386. True
  387. sage: Reducible_large_resp_sign_trop(f3,g,R2,w2)
  388. False
  389. """
  390. R0 = f[2].parent()
  391. flm = LM_trop(f[2],R,w)
  392. glm = LM_trop(g[2],R,w)
  393. nvar = len(f[2].parent().gens())
  394. trop_weight_sign = sum([f[1].degrees()[vv]*w[vv] for vv in range(nvar)])
  395. if glm.parent().monomial_divides(glm,flm):
  396. test_sign = R0.monomial_quotient(flm,glm)*g[1]
  397.  
  398. bool_comp_sign = (test_sign <= f[1])
  399. return (f[0]> g[0] and (test_sign.degree() <=f[1].degree())) or (f[0] == g[0] and (test_sign.degree() <f[1].degree())) or (f[0] == g[0] and (test_sign.degree() ==f[1].degree()) and bool_comp_sign)
  400. return False
  401.  
  402.  
  403. def Reducible_large_resp_sign_family_trop(f,G,R,w):
  404. r"""
  405.  
  406. Inside the F5 algorithm, f and the elements of the list G are of the form [i,x^alpha, poly].
  407. Tests whether G can reduce f while respecting their signatures
  408.  
  409. INPUT:
  410.  
  411. - ``f`` -- a polynomial with signature
  412. - ``G`` -- a list of polynomials with signature
  413. - ``R`` -- a polynomial ring of dimension
  414. - ``w`` -- a list of n rational numbers
  415.  
  416. OUTPUT:
  417.  
  418. a boolean, testing whether G can reduce f while respecting their signatures
  419.  
  420. EXAMPLES::
  421.  
  422. sage: S.<x,y,z>=QQ[]
  423. sage: R1.<x1,y1,z1> = Qp(2)[]
  424. sage: R2 = PolynomialRing(Qp(3),['x2','y2','z2'] , order = 'lex')
  425. sage: w1 = [0,0,0]
  426. sage: w2 = [0,0,-1]
  427. sage: from sage.rings.polynomial.padics.F5Trop_doctest import *
  428. sage: G = [[0,S(1),x^2], [1,S(1),x*y], [2,S(1),y*z], [3,S(1),z^2]]
  429. sage: f = [1,x^2,x^3*y]
  430. sage: f2 = [0,x^2,x^3]
  431. sage: f3 = [0,S(1),x^3]
  432. sage: Reducible_large_resp_sign_family_trop(f,G,R1,w1)
  433. True
  434. sage: Reducible_large_resp_sign_family_trop(f2,G,R1,w1)
  435. True
  436. sage: Reducible_large_resp_sign_family_trop(f3,G,R2,w2)
  437. False
  438. """
  439. for g in G:
  440. if Reducible_large_resp_sign_trop(f,g,R,w):
  441. return True
  442. return False
  443.  
  444. def Rewritten(R,G,g, R1,w):
  445. r"""
  446.  
  447. Find the latest element of G such it has a multiple with same signature as g.
  448.  
  449.  
  450. INPUT:
  451.  
  452. - ``R`` -- a polynomial ring
  453. - ``G`` -- a list of polynomials with signature, ordered by signature
  454. - ``g`` -- a polynomial with signature
  455. - ``R1`` -- a polynomial ring
  456. - ``G`` -- a list of polynomials with signature
  457.  
  458.  
  459. OUTPUT:
  460.  
  461. A multiple of an element of G with same signature and leading monomial
  462. as g. It is obtained with the latest element of g possible
  463.  
  464. EXAMPLES::
  465.  
  466. sage: S.<x,y,z>=QQ[]
  467. sage: R1.<x1,y1,z1> = Qp(2)[]
  468. sage: w1 = [0,0,0]
  469. sage: from sage.rings.polynomial.padics.F5Trop_doctest import *
  470. sage: G = [[0,y,x^2+2*z^2], [1,x,x*y+y^2], [2,1,y*z+z^2], [3,1,z^2]]
  471. sage: g = [2,x,2*x*y*z+x*z^2+y^2*z]
  472. sage: Rewritten(S,G,g,R1,w1)
  473. [2, x, x*y*z + x*z^2]
  474.  
  475. """
  476. mon_sign = LM_trop(R(g[1]),R1,w)
  477. newg=g
  478. for gg in G:
  479. if gg[0]==g[0] and R.monomial_divides(LM_trop(R(gg[1]),R1,w),mon_sign):
  480. glm = LM_trop(R(gg[1]),R1,w)
  481. mult = R.monomial_quotient(mon_sign,glm)
  482. newg = [gg[0],mult*gg[1],mult*gg[2]]
  483. return newg
  484.  
  485. def SymbolicPreprocessingTrop(R,w,G, d, paires,list_degrees):
  486. r"""
  487.  
  488. Constructs the Macaulay matrix from the pairs in paires,
  489. while using the F5 criterion
  490.  
  491.  
  492. INPUT:
  493.  
  494. - ``R`` -- a polynomial ring
  495. - ``w`` -- a list of n rational numbers
  496. - ``G`` -- a list of polynomials with signature
  497. - ``d`` -- an integer, representing the (sugar) degree of the processed polynomials
  498. - ``paires`` -- a list of pairs of polynomials with signature
  499. - ``list_degrees`` -- a list of integers, representing the degrees of the polynomials in F
  500.  
  501.  
  502. OUTPUT:
  503.  
  504. A Macaulay matrix, as a 3-tuple of a matrix built from polynomials,
  505. a list of monomials corresponding to its columns and
  506. a list of signatures for the rows of this matrix
  507.  
  508. EXAMPLES::
  509.  
  510. sage: S.<x,y,z>=QQ[]
  511. sage: R1.<x1,y1,z1> = Qp(2)[]
  512. sage: w1 = [0,0,0]
  513. sage: from sage.rings.polynomial.padics.F5Trop_doctest import*
  514. sage: G = [[0,y,x^2+2*z^2], [1,x,x*y+y^2], [2,S(1),y*z+z^2], [3,S(1),z^2]]
  515. sage: paires = []
  516. sage: paires.append([3,1,x,[1,x,x*y+y^2],y, [0,S(1),x^2+2*z^2]])
  517. sage: paires.append([3,3,y,[3,S(1),x*z+3*y*z+z^2],z, [3,S(1),x*y+2*y*z+z^2]])
  518. sage: d = 3
  519. sage: Mac = SymbolicPreprocessingTrop(R1,w1, G, d, paires,[1,1,2,2])
  520. sage: Mac[0][0]
  521. [1 1]
  522. [1 0]
  523. [0 1]
  524. sage: Mac[0][1]
  525. [[2, z], [3, z], [3, y]]
  526. sage: Mac[1]
  527. [z^3, y*z^2]
  528.  
  529. """
  530. R0 = G[0][2].parent()
  531.  
  532. paires_to_process = []
  533. paires_to_remove = []
  534. list_already_erased = []
  535. list_already_added = []
  536. # Might be faster with a list like list_already_passing_F5
  537.  
  538. for paire in paires:
  539. i_acceptable = True
  540. j_acceptable = True
  541. if paire[0] == d :
  542. paires_to_remove.append(paire)
  543. signj = [paire[3][0],paire[3][1]*paire[2]]
  544. signi = [paire[5][0],paire[5][1]*paire[4]]
  545. if not(signj in list_already_erased) and not(signi in list_already_erased):
  546. if F5_criterion_trop_affine(signj,G,R,w,list_degrees):
  547. j_acceptable = False
  548. list_already_erased.append(signj)
  549.  
  550. if F5_criterion_trop_affine(signi,G,R,w,list_degrees):
  551. i_acceptable = False
  552. list_already_erased.append(signi)
  553.  
  554. if i_acceptable and j_acceptable:
  555. paires_to_process.append(paire)
  556.  
  557.  
  558.  
  559. #Erasing the pairs that we have just considered
  560. for paire in paires_to_remove:
  561. paires.remove(paire)
  562.  
  563. #Producing the list of polynomials to be processed
  564. # We produce a set of the leading monomials
  565.  
  566. poly_to_do = []
  567. list_poly = []
  568. set_done_monom = set()
  569. for paire in paires_to_process:
  570. signj_polyj = [paire[3][0],paire[3][1]*paire[2],paire[2]*paire[3][2]]
  571. signj = [paire[3][0],paire[3][1]*paire[2]]
  572. if poly_to_do.count(signj_polyj)<1 and list_already_added.count(signj)<1:
  573. poly_to_do.append(signj_polyj)
  574. list_poly.append(signj_polyj)
  575. list_already_added.append(signj)
  576.  
  577. #set_done_monom.add(LM_trop(signj_polyj[2],R,w))
  578.  
  579. signi_polyi = [paire[5][0],paire[5][1]*paire[4], paire[4]*paire[5][2]]
  580. signi = [paire[5][0],paire[5][1]*paire[4]]
  581. if poly_to_do.count(signi_polyi)<1 and list_already_added.count(signi)<1:
  582. poly_to_do.append(signi_polyi)
  583. list_poly.append(signi_polyi)
  584. list_already_added.append(signi)
  585.  
  586.  
  587.  
  588. #Here, we add the initial polynomials that are of degree d to the matrix.
  589. #The goal is to prevent redundancy that might happen when an initial polynomial is far from being reduced
  590. for g in G:
  591. if g[2].degree() == d and g[1] == 1 and poly_to_do.count(g)<1 and list_already_added.count([g[0],g[1]])<1 and not([g[0],g[1]] in list_already_erased):
  592. poly_to_do.append(g)
  593. list_poly.append(g)
  594. list_already_added.append([g[0],g[1]])
  595.  
  596.  
  597.  
  598. poly_to_do.sort()
  599. list_poly.sort()
  600.  
  601. #Rewritting
  602.  
  603. len_poly_to_do = len(poly_to_do)
  604. for ss in range(len_poly_to_do):
  605. poly_to_do[ss] = Rewritten(R0,G, poly_to_do[ss],R,w)
  606. list_poly[ss] = poly_to_do[ss]
  607.  
  608. monom = set([])
  609.  
  610.  
  611. while len(poly_to_do)>0:
  612. sign_et_poly = poly_to_do.pop(0)
  613. monom = monom.union(set(sign_et_poly[2].monomials()))
  614. set_monom = set(sign_et_poly[2].monomials())-set_done_monom
  615. set_monom = list(set_monom)
  616. set_monom.sort()
  617.  
  618. while len(set_monom)>0:
  619. mon = set_monom.pop()
  620.  
  621. # We do take for mon the biggest monomial available here
  622. # We look for a g in G and x^alpha such that mon ==x^alpha* g.lm() with smallest signature
  623. # In other word, we look for a "reducer." Yet for signature reason, this reducer could be more reduced than reducing
  624. smallest_g = 0
  625. smallest_sign = []
  626. for g in G:
  627. if mon.parent().monomial_divides(LM_trop(g[2],R,w),mon):
  628. glm = LM_trop(g[2],R,w)
  629. mult = mon.parent().monomial_quotient(mon,glm)
  630. sign_to_test = [g[0],mult*g[1]]
  631. if len(smallest_sign)==0 and list_already_added.count(sign_to_test)<1:
  632. if not(F5_criterion_trop_affine(sign_to_test,G,R,w,list_degrees)):
  633. smallest_g = mult*g[2]
  634. smallest_sign = [g[0],mult*g[1]]
  635. else:
  636. if [g[0],mult*g[1]] <= smallest_sign and list_already_added.count([g[0],mult*g[1]])<1:
  637. if not(F5_criterion_trop_affine(sign_to_test,G,R,w,list_degrees)):
  638. smallest_g = mult*g[2]
  639. smallest_sign = [g[0],mult*g[1]]
  640.  
  641. if len(smallest_sign)>0:
  642. list_already_added.append(smallest_sign)
  643. new_poly = smallest_sign+[smallest_g]
  644. poly_to_do.append(new_poly)
  645. poly_to_do.sort()
  646. list_poly.append(new_poly)
  647. list_poly.sort()
  648. set_done_monom.add(mon)
  649.  
  650. monom = list(monom)
  651. monom.sort()
  652.  
  653.  
  654.  
  655. if len(list_poly)>0:
  656. Mac = Make_Macaulay_Matrix_sparse_dict(list_poly, monom)
  657. else:
  658. Mac = [Matrix(R.base_ring(),0,0,[]),[]]
  659.  
  660. return Mac, monom, paires
  661.  
  662. def list_indexes_trop(Mat, monom, R,w):
  663. r"""
  664.  
  665. Computes the list of the indexes of the leading coefficients
  666. of the rows of the Macaulay matrix Mat, given with
  667. the list of monomials monom indexing its columns
  668. and R and w defining a tropical term order on the
  669. polynomial ring in which the monomials of monom live.
  670.  
  671.  
  672. INPUT:
  673.  
  674. - ``Mat`` -- a matrix
  675. - ``monom`` -- a list of monomials
  676. - ``R`` -- a polynomial ring of dimension n
  677. - ``w`` -- a list of n rational numbers
  678.  
  679. OUTPUT:
  680.  
  681. a list of integers
  682.  
  683. EXAMPLES::
  684.  
  685. sage: S.<x,y,z,t>=QQ[]
  686. sage: from sage.rings.polynomial.padics.F5Trop_doctest import *
  687. sage: R1.<x1,y1,z1,t1> = Qp(2)[]
  688. sage: w1 = [0,0,0,0]
  689. sage: mat = Matrix(QQ,4,4,[0,1,2,0,0,0,0,3,2,0,0,0,0,0,1,1])
  690. sage: monom = [x,y,z,t]
  691. sage: list_indexes_trop(mat, monom, R1, w1)
  692. [1, 3, 0, 2]
  693. sage: w2 = [2,9,1,2]
  694. sage: list_indexes_trop(mat, monom, R1, w2)
  695. [2, 3, 0, 2]
  696. """
  697. n = Mat.nrows()
  698. m = Mat.ncols()
  699. listindexes = []
  700. T = R.base_ring()
  701. K=T
  702. nvar = len(w)
  703.  
  704. def to_R(mon):
  705. d = PolyDict({ETuple(mon.exponents()[0]): K.one()}, force_etuples=False)
  706. return MPolynomial_polydict(R, d)
  707.  
  708. for i in range(n):
  709. #Looking for the LT for the tropical order
  710. max_term = None
  711. max_col = None
  712. entry_non_zero = None
  713. max_deg = None
  714. for j in range(m):
  715. if max_term is None and Mat[i,j] != 0:
  716. max_term = Mat[i,j]
  717. max_col = j
  718. max_monomial = monom[j]
  719. max_deg = monom[j].degree()
  720. max_trop_weight = T(Mat[i,j]).valuation()+sum([monom[j].degrees()[uu]*w[uu] for uu in range(nvar)])
  721. elif Mat[i,j] != 0:
  722. deg = monom[j].degree()
  723. term_trop_weight = T(Mat[i,j]).valuation()+sum([monom[j].degrees()[uu]*w[uu] for uu in range(nvar)])
  724. if deg>max_deg or (deg== max_deg and term_trop_weight < max_trop_weight) or (deg== max_deg and term_trop_weight == max_trop_weight and to_R(monom[j])> to_R(max_monomial)):
  725. max_term = Mat[i,j]
  726. max_col = j
  727. max_trop_weight = term_trop_weight
  728. max_deg = deg
  729. max_monomial = monom[j]
  730. if max_term is None:
  731. listindexes.append(m)
  732. else:
  733. listindexes.append(max_col)
  734. return listindexes
  735.  
  736. def column_transposition(mat,list_monom,j1,j2):
  737. """
  738. performs the transposition of the columns j1 and j2 of the Macaulay matrix given by mat and list_monom.
  739.  
  740. EXAMPLES::
  741.  
  742. sage: from sage.rings.polynomial.padics.F5Trop_doctest import column_transposition
  743. sage: T.<x,y,z> = QQ[]
  744. sage: Mac = Matrix(QQ,1,3,[1,2,3])
  745. sage: list_monom = [x,y,z]
  746. sage: column_transposition(Mac,list_monom,1,2)
  747. sage: Mac
  748. [1 3 2]
  749. sage: list_monom
  750. [x, z, y]
  751.  
  752. """
  753. if j1 != j2:
  754. nrows = mat.nrows()
  755. ncols = mat.ncols
  756. temp = list_monom[j1]
  757. list_monom[j1] = list_monom[j2]
  758. list_monom[j2] = temp
  759. mat.swap_columns(j1,j2)
  760.  
  761.  
  762. def Trop_RowEchelonMac_and_listpivots(Mac,list_monom0, R, tropical_weight):
  763. r"""
  764.  
  765. Computes the row-echelon form of the Macaulay matrix Mac,
  766. with no choice of pivot: always the LM of each row is used as pivot
  767.  
  768.  
  769. INPUT:
  770.  
  771. - ``Mac`` -- a matrix
  772. - ``list_monom`` -- a list of the monomials of degree d in R
  773. - ``R`` -- a polynomial ring
  774. - ``tropical_weight`` -- a tropical weight, i.e. an n-tuple of real numbers
  775.  
  776. OUTPUT:
  777.  
  778. a list consisting of Mactilde, a matrix under row-echelon form
  779. (up to permutation of its rows), a list of the
  780. indexes of the pivots used for the computation of
  781. Mactilde (number of rows plus one if there is no
  782. pivot on a row) and a list
  783.  
  784. EXAMPLES::
  785.  
  786. sage: from sage.rings.polynomial.padics.F5Trop_doctest import *
  787. sage: S.<x,y,z,t> = QQ[]
  788. sage: R.<xx,yy,zz,tt> = Qp(2,50)[]
  789. sage: w = [-1,0,2,-3]
  790. sage: list_monom = [x,y,z,t]
  791. sage: mat = Matrix(QQ,4,4,[0,1,2,0,-1,-2,0,3,2,0,0,0,2,1,3,1])
  792. sage: u = Trop_RowEchelonMac_and_listpivots(mat, list_monom,R,w)
  793. sage: u[0]
  794. [ 1 0 0 2]
  795. [ 0 1 -1/3 4/3]
  796. [ 0 0 1 0]
  797. [ 0 0 0 1]
  798. sage: u[1]
  799. [y, t, x, z]
  800.  
  801. """
  802. n = Mac.nrows()
  803. m = Mac.ncols()
  804. Mactilde = copy(Mac)
  805. list_monom = [u for u in list_monom0]
  806. listpivots = []
  807. T = R.base_ring()
  808. K=T
  809. def to_R(mon):
  810. d = PolyDict({ETuple(mon.exponents()[0]): K.one()}, force_etuples=False)
  811. return MPolynomial_polydict(R, d)
  812.  
  813. nvar = len(tropical_weight)
  814. index_column = 0
  815. for i in range(n):
  816. #Looking for the LT for the tropical order
  817. max_term = None
  818. max_col = None
  819. entry_non_zero = None
  820. max_deg = None
  821. for j0 in range(index_column,m):
  822. if max_term is None and Mactilde[i,j0] != 0:
  823. max_term = Mactilde[i,j0]
  824. max_col = j0
  825. max_monomial = list_monom[j0]
  826. max_deg = list_monom[j0].degree()
  827. max_trop_weight = T(Mactilde[i,j0]).valuation()+sum([list_monom[j0].degrees()[uu]*tropical_weight[uu] for uu in range(nvar)])
  828. elif Mactilde[i,j0] != 0:
  829. term_trop_weight = T(Mactilde[i,j0]).valuation()+sum([list_monom[j0].degrees()[uu]*tropical_weight[uu] for uu in range(nvar)])
  830. deg = list_monom[j0].degree()
  831. if deg>max_deg or (deg== max_deg and term_trop_weight < max_trop_weight) or (deg== max_deg and term_trop_weight == max_trop_weight and to_R(list_monom[j0])> to_R(max_monomial)):
  832. max_term = Mactilde[i,j0]
  833. max_col = j0
  834. max_trop_weight = term_trop_weight
  835. max_deg = deg
  836. max_monomial = list_monom[j0]
  837. if max_term is None:
  838. listpivots.append(m)
  839. else:
  840. listpivots.append(index_column)
  841. column_transposition(Mactilde, list_monom, index_column, max_col)
  842. Mactilde.rescale_row(i, Mactilde[i,index_column]**(-1), start_col=index_column)
  843. for l in range(i+1,n):
  844. if Mactilde[l,index_column] != 0:
  845. Mactilde.add_multiple_of_row(l, i, -Mactilde[l,index_column] , start_col=index_column)
  846. index_column +=1
  847. return Mactilde, list_monom, listpivots
  848.  
  849.  
  850. def Eliminate_unreduced(G,F,d):
  851. r"""
  852.  
  853. In the list of polynomials with signature G, remove the polynomials
  854. of sugar degree d that are both in G and F, leading to a redundancy of signature
  855.  
  856.  
  857. INPUT:
  858.  
  859. - ``G`` -- a list of polynomials with signature, ordered by signature
  860. - ``F`` -- a list of polynomials with signature, ordered by signature
  861. - ``d`` -- an integer, to tell the degree of the polynomials that are processed
  862.  
  863. OUTPUT:
  864.  
  865. - ``newG`` -- an updated list of polynomials with signature
  866.  
  867. EXAMPLES::
  868.  
  869. sage: from sage.rings.polynomial.padics.F5Trop_doctest import *
  870. sage: S.<x,y,z>=QQ[]
  871. sage: F = [[0,1,x],[1,1,x^2+2*x*y+y^2-z^2]]
  872. sage: G = [[0,1,x],[1,1,x^2+2*x*y+y^2-z^2],[1,1,y^2-z^2]]
  873. sage: G=Eliminate_unreduced(G,F,2)
  874. sage: G
  875. [[0, 1, x], [1, 1, y^2 - z^2]]
  876. """
  877.  
  878. newGd = []
  879. newG = []
  880. newGsignd = []
  881. for g in G:
  882. newG.append(g)
  883. if g[2].degree() == d:
  884. newGd.append(g)
  885. newGsignd.append([g[0],g[1]])
  886. tempF = []
  887. tempFsign = []
  888. for f in F:
  889. if f[2].degree() == d:
  890. tempF.append(f)
  891. tempFsign.append([f[0],f[1]])
  892. dl = len(tempF)
  893. for i in range(dl):
  894. if newGsignd.count(tempFsign[i])>1:
  895. if newGd.count(tempF[i])>0:
  896. newG.remove(tempF[i])
  897. return newG
  898.  
  899. def UpdateGpaires_trop(Mac, monom, Mactilde, lmonom_mactilde, listpivots, G, paires, R1, w, redundancy_test, F,d, list_degrees):
  900. r"""
  901.  
  902. In the F5 algorithm, updates the current list of polynomials with
  903. signatures G, according to the Macaulay matrix Mac,
  904. its row-echelon form (with no choice of pivot) Mactilde,
  905. and previous computation.
  906. paires, the list of the remaining S-pairs to proceed is also
  907. updated.
  908. G and paires are updated in place.
  909.  
  910.  
  911. INPUT:
  912.  
  913. - ``Mac`` -- a Macaulay matrix: a list of a matrix and a list of the signatures corresponding to its rows
  914. - ``monom`` -- the list of the monomials corresponding to the columns of Mac
  915. - ``Mactilde`` -- a Macaulay matrix (just the matrix), tropical row-echelon form of Mac
  916. - ``lmonom_mactilde`` -- the list of the monomials corresponding to the columns of Mac
  917. - ``listpivots`` -- the list of the pivots used when computing the tropical row-echelon form Mactilde of Mac
  918. - ``G`` -- a list of polynomials with signature
  919. - ``paires`` -- a list of S-pairs of polynomials with signature
  920. - ``R1`` -- a polynomial ring
  921. - ``w`` -- a tropical weight, i.e. an n-tuple of real numbers
  922. - ``redundancy_test`` -- a boolean to check whether a redundancy test is required
  923. - ``F`` -- a list of polynomials with signature, ordered by signature
  924. - ``d`` -- an integer, to tell the degree of the polynomials that are processed
  925. - ``list_degrees`` -- a list of integers, representing the degrees of the polynomials in F
  926.  
  927. OUTPUT:
  928.  
  929. - ``G`` -- an updated list of polynomials with signature
  930. - ``paires`` -- an updated list of S-pairs of polynomials with signature
  931.  
  932. EXAMPLES::
  933.  
  934. sage: from sage.rings.polynomial.padics.F5Trop_doctest import *
  935. sage: S.<x,y,z>=QQ[]
  936. sage: R.<x1,y1,z1>=Qp(2)[]
  937. sage: w = [0,0,0]
  938. sage: mat = Matrix(QQ,2,6,[2,1,0,0,0,4,0,1,1,0,0,-1])
  939. sage: monom = [x^2,x*y,y^2,x*z,y*z,z^2]
  940. sage: Mac = [mat,[[0,S(1)],[1,S(1)]]]
  941. sage: mat2 = Matrix(QQ,2,6,[1,0,2,0,0,4,0,1,-2,0,0,-9])
  942. sage: lmonom_mat2 = [x*y,y^2, x^2,x*z,y*z,z^2]
  943. sage: listpivots = [0,1]
  944. sage: G = [[0,1,2*x^2+x*y+4*z^2],[1,1,x*y+y^2-z^2]]
  945. sage: F=G
  946. sage: paires = []
  947. sage: d = 2
  948. sage: G, paires = UpdateGpaires_trop(Mac, monom, mat2, lmonom_mat2, listpivots, G, paires, R,w, false, F,d,[2,2])
  949. sage: G
  950. [[0, 1, 2*x^2 + x*y + 4*z^2],
  951. [1, 1, x*y + y^2 - z^2],
  952. [1, 1, -2*x^2 + y^2 - 9*z^2]]
  953. sage: paires
  954. [[3, 1, x, [1, 1, -2*x^2 + y^2 - 9*z^2], y, [0, 1, 2*x^2 + x*y + 4*z^2]],
  955. [3, 1, x, [1, 1, -2*x^2 + y^2 - 9*z^2], y, [1, 1, x*y + y^2 - z^2]]]
  956.  
  957. """
  958. n = Mac[0].nrows()
  959. m = Mac[0].ncols()
  960.  
  961. if n == 0:
  962. return G, paires
  963.  
  964. listindexes = list_indexes_trop(Mac[0], monom,R1,w)
  965. R = G[0][2].parent()
  966.  
  967. newG = []
  968. for i in range(n):
  969. if listpivots[i]<m and monom[listindexes[i]]<>lmonom_mactilde[listpivots[i]]:
  970. if not(Reducible_large_resp_sign_family_trop([Mac[1][i][0],Mac[1][i][1],lmonom_mactilde[listpivots[i]]],G,R1,w)):
  971. newG.append(Reconstruction(Mactilde,i,lmonom_mactilde,[Mac[1][i][0],Mac[1][i][1]]))
  972.  
  973. G2 = G+newG
  974. G2.sort()
  975.  
  976. #Redundancy elimination
  977. G3= []
  978. for g in G2:
  979. if G3.count(g)==0:
  980. G3.append(g)
  981. G2=G3
  982.  
  983. if redundancy_test:
  984. paires = []
  985. G2=Eliminate_unreduced(G2,F,d)
  986. G = [aa for aa in G if G2.count(aa)>0]
  987. newG = [aa for aa in newG if G2.count(aa)>0]
  988.  
  989. long = len(G2)
  990. for j1 in range(long):
  991. for j2 in range(j1+1,long):
  992. g = G2[j1]
  993. g2 = G2[j2]
  994. i1 = g[0]
  995. i2 = g2[0]
  996. gtemp1 = g
  997. gtemp2 = g2
  998. if i2< i1:
  999. h = gtemp1
  1000. gtemp1 = gtemp2
  1001. gtemp2 = h
  1002. u = R.monomial_lcm(LM_trop(gtemp1[2],R1,w),LM_trop(gtemp2[2],R1,w))
  1003. quot = R.monomial_quotient(u,LM_trop(gtemp2[2],R1,w))
  1004. sugar = (quot*gtemp2[1]).degree()+list_degrees[gtemp2[0]]
  1005. if sugar > d:
  1006. paires.append([sugar,gtemp2[0],R.monomial_quotient(u,LM_trop(gtemp2[2],R1,w)),gtemp2,R.monomial_quotient(u,LM_trop(gtemp1[2],R1,w)), gtemp1])
  1007.  
  1008. else:
  1009.  
  1010. for g in G:
  1011. for g2 in newG:
  1012. i1 = g[0]
  1013. i2 = g2[0]
  1014. gtemp1 = g
  1015. gtemp2 = g2
  1016. if i2< i1:
  1017. h = gtemp1
  1018. gtemp1 = gtemp2
  1019. gtemp2 = h
  1020. u = R.monomial_lcm(LM_trop(gtemp1[2],R1,w),LM_trop(gtemp2[2],R1,w))
  1021. quot = R.monomial_quotient(u,LM_trop(gtemp2[2],R1,w))
  1022. sugar = (quot*gtemp2[1]).degree()+list_degrees[gtemp2[0]]
  1023. if sugar > d:
  1024. paires.append([sugar,gtemp2[0],R.monomial_quotient(u,LM_trop(gtemp2[2],R1,w)),gtemp2,R.monomial_quotient(u,LM_trop(gtemp1[2],R1,w)), gtemp1])
  1025.  
  1026. long = len(newG)
  1027. for j1 in range(long):
  1028. for j2 in range(j1+1,long):
  1029. g = newG[j1]
  1030. g2 = newG[j2]
  1031. i1 = g[0]
  1032. i2 = g2[0]
  1033. gtemp1 = g
  1034. gtemp2 = g2
  1035. if i2< i1:
  1036. h = gtemp1
  1037. gtemp1 = gtemp2
  1038. gtemp2 = h
  1039. u = R.monomial_lcm(LM_trop(gtemp1[2],R1,w),LM_trop(gtemp2[2],R1,w))
  1040. quot = R.monomial_quotient(u,LM_trop(gtemp2[2],R1,w))
  1041. sugar = (quot*gtemp2[1]).degree()+list_degrees[gtemp2[0]]
  1042. if sugar > d:
  1043. paires.append([sugar,gtemp2[0],R.monomial_quotient(u,LM_trop(gtemp2[2],R1,w)),gtemp2,R.monomial_quotient(u,LM_trop(gtemp1[2],R1,w)), gtemp1])
  1044.  
  1045.  
  1046. if len(newG)>0 or redundancy_test:
  1047. paires.sort()
  1048.  
  1049. return G2, paires
  1050.  
  1051.  
  1052.  
  1053. def tentative_F5Trop(list1,R1,w):
  1054. r"""
  1055. Computes a tropical Gröbner bases of the ideal generated by the list of polynomials list1.
  1056. The tropical term order is defined by the valuation on the base ring of R1, its ambiant
  1057. monomial order, and the tropical weight w.
  1058.  
  1059.  
  1060. INPUT:
  1061.  
  1062. - ``list1`` -- a list of polynomials
  1063. - ``R1`` -- a polynomial ring
  1064. - ``w`` -- a tropical weight, i.e. an n-tuple of real numbers
  1065.  
  1066. OUTPUT:
  1067.  
  1068. A Gröbner basis of the ideal generated by list1. For now, it is not reduced
  1069.  
  1070.  
  1071. EXAMPLES::
  1072.  
  1073. sage: from sage.rings.polynomial.padics.F5Trop_doctest import *
  1074. sage: S.<x,y,z>=QQ[]
  1075. sage: R.<x1,y1,z1>=Qp(2)[]
  1076. sage: w = [0,0,0]
  1077. sage: list1 = [x+y+z, x*y+x*z+y*z,x*y*z]
  1078. sage: G = tentative_F5Trop(list1,R,w)
  1079. sage: G
  1080. [x + y + z, y^2 + y*z + z^2, z^3]
  1081.  
  1082. """
  1083. #Initialisation
  1084. R = list1[0].parent()
  1085. list_degrees = [aa.degree() for aa in list1]
  1086. list = copy(list1)
  1087. s = len(list)
  1088. listsign = []
  1089. for i in range(s):
  1090. listsign.append([i,R(1),list[i]])
  1091. G = listsign
  1092. l = len(G)
  1093. paires = []
  1094. for i in range(s):
  1095. for j in range(i+1,l):
  1096.  
  1097. u = R.monomial_lcm(LM_trop(list[i],R1,w),LM_trop(list[j],R1,w))
  1098. quot = R.monomial_quotient(u,LM_trop(list[j],R1,w))
  1099. sugar = (quot*R(1)).degree()+list_degrees[j]
  1100. paires.append([sugar,j,R.monomial_quotient(u,LM_trop(list[j],R1,w)),listsign[j],R.monomial_quotient(u,LM_trop(list[i],R1,w)), listsign[i]])
  1101. #watch out: j>i
  1102.  
  1103. paires.sort()
  1104. #Main part
  1105. d = paires[0][0]
  1106. while paires != [] :
  1107. # print "d ="
  1108. # print(d)
  1109. # print "paires restantes"
  1110. # print len(paires)
  1111. Mac, monom, paires = SymbolicPreprocessingTrop(R1,w,G, d,paires, list_degrees)
  1112. Mactilde, list_monom_mactilde, listpivots = Trop_RowEchelonMac_and_listpivots(Mac[0],monom, R1, w)
  1113. # print Mactilde.nrows()
  1114. # print monom
  1115. # print Mac[0].str()
  1116. # print Mac[1]
  1117. # print "resultat de la reduction"
  1118. # print list_monom_mactilde
  1119. # print Mactilde.str()
  1120. redundancy_check = (list_degrees.count(d)>0)
  1121. G, paires = UpdateGpaires_trop(Mac, monom, Mactilde,list_monom_mactilde, listpivots, G, paires, R1,w,redundancy_check, listsign,d, list_degrees)
  1122. # print "G"
  1123. # print G
  1124. # print "les lm"
  1125. # print [[g[0],g[1], LM_trop(g[2],R1,w)] for g in G]
  1126. # print "les paires"
  1127. # print paires[0]
  1128. d=d+1
  1129. # print "G final"
  1130. # print G
  1131. # print "les lm"
  1132. # print [[g[0],g[1], LM_trop(g[2],R1,w)] for g in G]
  1133. #print "G final"
  1134. #print G
  1135. G = [g[2] for g in G]
  1136. return G
  1137.  
  1138.  
  1139. ##################################################
  1140. # Classical reduction for computing a reduced GB #
  1141. ##################################################
  1142.  
  1143. def minimalization(G,R1,w):
  1144. r"""
  1145. Computes a minimal tropical Gröbner basis of the ideal generated by the list of polynomials in G.
  1146. The tropical term order is defined by the valuation on the base ring of R1, its ambiant
  1147. monomial order, and the tropical weight w.
  1148. G is a tropical Gröbner basis (for the same term order).
  1149.  
  1150.  
  1151. INPUT:
  1152.  
  1153. - ``G`` -- a list of polynomials, which is a tropical Gröbner basis
  1154. - ``R1`` -- a polynomial ring
  1155. - ``w`` -- a tropical weight, i.e. an n-tuple of real numbers
  1156.  
  1157. OUTPUT:
  1158.  
  1159. A minimal tropical Gröbner basis of the ideal generated by G. For now, it is not reduced
  1160.  
  1161.  
  1162. EXAMPLES::
  1163.  
  1164. sage: from sage.rings.polynomial.padics.F5Trop_doctest import *
  1165. sage: S.<x,y,z>=QQ[]
  1166. sage: R.<x1,y1,z1>=Qp(2)[]
  1167. sage: w = [0,0,0]
  1168. sage: list1 = [x+y+z, x*y+x*z+y*z,x*y*z]
  1169. sage: G = tentative_F5Trop(list1,R,w)
  1170. sage: minimalization(G,R,w)
  1171. [x + y + z, y^2 + y*z + z^2, z^3]
  1172.  
  1173. """
  1174. P = G[0].parent()
  1175. list_minimal = []
  1176. list_lm_minimal = []
  1177. listlm = []
  1178. for g in G:
  1179. listlm.append(LM_trop(g,R1,w))
  1180. i = 0
  1181. r = len(G)
  1182. while i < len(listlm):
  1183. flag = True
  1184. a = listlm[i]
  1185. if list_lm_minimal.count(a)>0:
  1186. flag = False
  1187. j = 0
  1188. while (j < r and flag):
  1189. b = listlm[j]
  1190. if (a != b) and P.monomial_divides(b,a):
  1191. flag = False
  1192. j = j + 1
  1193. if flag:
  1194. list_minimal.append(G[i])
  1195. list_lm_minimal.append(a)
  1196. i = i + 1
  1197. return list_minimal
  1198.  
  1199. ###############################################
  1200. # Reduction of a Trop GB #
  1201. ###############################################
  1202.  
  1203. def symb_preprocess_for_reduction(G,R1,w):
  1204. r"""
  1205. Computes a list of 3-tuples (integer, integer, Macaulay matrix) , starting from G a minimal tropical GB, in order to do the inter-reduction of G
  1206.  
  1207.  
  1208. INPUT:
  1209.  
  1210. - ``G`` -- a list of polynomials, which is a minimal tropical Gröbner basis
  1211. - ``R1`` -- a polynomial ring
  1212. - ``w`` -- a tropical weight, i.e. an n-tuple of real numbers
  1213.  
  1214. OUTPUT:
  1215.  
  1216. A list of couples 3-tuples (integer, integer, Macaulay matrix)
  1217. The first one is the degree of the Macaulay matrix, the second is the number of rows coming from
  1218. polynomials of G in the matrix (placed at the top)
  1219. a list indexed by d of the list of monomials of degree d, for all d up to the maximal required
  1220.  
  1221.  
  1222. EXAMPLES::
  1223.  
  1224. sage: from sage.rings.polynomial.padics.F5Trop_doctest import *
  1225. sage: S.<x,y,z>=QQ[]
  1226. sage: R.<x1,y1,z1>=Qp(2)[]
  1227. sage: w = [0,0,0]
  1228. sage: list1 = [z^2, x^3 + 2*x^2*y + y^3, y^3 + y^2*z + y*z^2, y^2*z + z^3]
  1229. sage: G = tentative_F5Trop(list1,R,w)
  1230. sage: G = minimalization(G,R,w)
  1231. sage: u = symb_preprocess_for_reduction(G,R,w)
  1232. sage: u[0]
  1233. [0 0 0 0 1 0 0]
  1234. [1 0 2 0 0 0 1]
  1235. [1 1 0 1 0 0 0]
  1236. [0 0 0 1 0 1 0]
  1237. [1 1 0 1 0 0 0]
  1238. [0 0 0 1 0 1 0]
  1239. [0 1 0 0 0 0 0]
  1240. [0 0 0 0 0 1 0]
  1241. sage: u[1]
  1242. [y^3, y*z^2, x^2*y, y^2*z, z^2, z^3, x^3]
  1243. """
  1244. lmonom = set([])
  1245. P = G[0].parent()
  1246.  
  1247. for g in G:
  1248. lmonom = lmonom.union(g.monomials())
  1249.  
  1250. list_to_do = []
  1251. list_red = []
  1252. set_done_monom = set([])
  1253. for g in G:
  1254. list_to_do.append(g)
  1255. set_done_monom.add(LM_trop(g,R1,w))
  1256. while len(list_to_do)>0:
  1257. poly = list_to_do.pop(0)
  1258. set_monom = set(poly.monomials())-set_done_monom
  1259. set_monom = list(set_monom)
  1260.  
  1261. while len(set_monom)>0:
  1262. mon = set_monom.pop()
  1263. for g in G:
  1264. glm = LM_trop(g,R1,w)
  1265. if mon.parent().monomial_divides(glm,mon):
  1266. mult = mon.parent().monomial_quotient(mon,glm)
  1267. list_red.append(mult*g)
  1268. list_to_do.append(mult*g)
  1269. lmonom = lmonom.union(set((mult*g).monomials()))
  1270. break
  1271. set_done_monom.add(mon)
  1272. listd=G+list_red
  1273. lmonom = list(lmonom)
  1274. mat = Make_Macaulay_Matrix_sparse_dict_no_sign(listd,lmonom)
  1275. return mat, lmonom
  1276.  
  1277. def tropical_complete_row_reduction(Mat, lmonom, R,w):
  1278. r"""
  1279. Computes a tropical complete row reduction of the Macaulay matrix Mat, with columns indexed by lmonom
  1280. and with term order defined by R1 and w
  1281.  
  1282.  
  1283. INPUT:
  1284.  
  1285. - ``Mat`` -- a matrix
  1286. - ``lmonom`` -- a list of the monomials of a given degree
  1287. - ``R`` -- a polynomial ring
  1288. - ``w`` -- a tropical weight, i.e. an n-tuple of real numbers
  1289.  
  1290. OUTPUT:
  1291.  
  1292. A reduced tropical Gröbner basis of the ideal generated by G.
  1293.  
  1294.  
  1295. EXAMPLES::
  1296.  
  1297. sage: from sage.rings.polynomial.padics.F5Trop_doctest import *
  1298. sage: S.<x,y,z>=QQ[]
  1299. sage: R.<x1,y1,z1>=Qp(2)[]
  1300. sage: w = [0,0,0]
  1301. sage: Mat = Matrix(QQ,5,10,[1,2,0,1,0,0,0,0,0,0,0,0,0,1,0,0,1,0,1,0,0,0,0,0,0,0,1,0,0,1,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,1])
  1302. sage: lmonom = [x^3,x^2*y,x*y^2,y^3,x^2*z,x*y*z,y^2*z, x*z^2,y*z^2,z^3]
  1303. sage: u = tropical_complete_row_reduction(Mat, lmonom, R,w)
  1304. sage: u[0]
  1305. [1 0 0 0 0 0 0 0 2 0]
  1306. [0 1 0 0 0 0 0 0 0 0]
  1307. [0 0 1 0 0 0 0 0 0 0]
  1308. [0 0 0 1 0 0 0 0 0 0]
  1309. [0 0 0 0 1 0 0 0 0 0]
  1310. sage: u[1]
  1311. [x^3, y^3, y^2*z, y*z^2, z^3, x*y*z, x*y^2, x*z^2, x^2*y, x^2*z]
  1312.  
  1313. """
  1314. n = Mat.nrows()
  1315. m = Mat.ncols()
  1316. Mactilde = copy(Mat)
  1317. list_monom = copy(lmonom)
  1318. T = R.base_ring()
  1319. K=T
  1320. def to_R(mon):
  1321. d = PolyDict({ETuple(mon.exponents()[0]): K.one()}, force_etuples=False)
  1322. return MPolynomial_polydict(R, d)
  1323. nvar = len(w)
  1324. index_column = 0
  1325. for i in range(n):
  1326. #Looking for the LT for the tropical order
  1327. max_term = None
  1328. max_col = None
  1329. entry_non_zero = None
  1330. max_deg = None
  1331. for j0 in range(index_column,m):
  1332. if max_term is None and Mactilde[i,j0] != 0:
  1333. max_term = Mactilde[i,j0]
  1334. max_col = j0
  1335. max_monomial = list_monom[j0]
  1336. max_deg = list_monom[j0].degree()
  1337. max_trop_weight = T(Mactilde[i,j0]).valuation()+sum([list_monom[j0].degrees()[uu]*w[uu] for uu in range(nvar)])
  1338. elif Mactilde[i,j0] != 0:
  1339. term_trop_weight = T(Mactilde[i,j0]).valuation()+sum([list_monom[j0].degrees()[uu]*w[uu] for uu in range(nvar)])
  1340. deg = list_monom[j0].degree()
  1341. if deg>max_deg or (deg== max_deg and term_trop_weight < max_trop_weight) or (deg== max_deg and term_trop_weight == max_trop_weight and to_R(list_monom[j0])> to_R(max_monomial)):
  1342. max_term = Mactilde[i,j0]
  1343. max_col = j0
  1344. max_trop_weight = term_trop_weight
  1345. max_deg = deg
  1346. max_monomial = list_monom[j0]
  1347. if max_term is not None:
  1348. column_transposition(Mactilde, list_monom, index_column, max_col)
  1349. Mactilde.rescale_row(i, Mactilde[i,index_column]**(-1), start_col=index_column)
  1350. for l in range(n):
  1351. if l <> i:
  1352. Mactilde.add_multiple_of_row(l, i, -Mactilde[l,index_column] , start_col=index_column)
  1353. index_column +=1
  1354. return Mactilde, list_monom
  1355.  
  1356.  
  1357.  
  1358. def reduction_trop_GB(G,R1,w):
  1359. r"""
  1360. Computes a reduced tropical Gröbner basis of the ideal generated by the list of polynomials in G.
  1361. The tropical term order is defined by the valuation on the base ring of R1, its ambiant
  1362. monomial order, and the tropical weight w.
  1363. G is a minimal tropical Gröbner basis (for the same term order).
  1364.  
  1365.  
  1366. INPUT:
  1367.  
  1368. - ``G`` -- a list of polynomials, which is a tropical Gröbner basis
  1369. - ``R1`` -- a polynomial ring
  1370. - ``w`` -- a tropical weight, i.e. an n-tuple of real numbers
  1371.  
  1372. OUTPUT:
  1373.  
  1374. A reduced tropical Gröbner basis of the ideal generated by G.
  1375.  
  1376.  
  1377. EXAMPLES::
  1378.  
  1379. sage: from sage.rings.polynomial.padics.F5Trop_doctest import *
  1380. sage: S.<x,y,z>=QQ[]
  1381. sage: R.<x1,y1,z1>=Qp(2)[]
  1382. sage: w = [0,0,0]
  1383. sage: list1 = [z^2, x^3 + 2*x^2*y + y^3, y^3 + y^2*z + y*z^2, y^2*z + z^3]
  1384. sage: G = tentative_F5Trop(list1,R,w)
  1385. sage: G = minimalization(G,R,w)
  1386. sage: G = reduction_trop_GB(G,R,w)
  1387. sage: G
  1388. [z^2, x^3 + 2*x^2*y, y^3, y^2*z]
  1389.  
  1390. """
  1391. mat, lmonom = symb_preprocess_for_reduction(G,R1,w)
  1392. P = G[0].parent()
  1393. G2 = []
  1394. nb = len(G)
  1395. Mactilde, lmonom2 = tropical_complete_row_reduction(mat, lmonom, R1,w)
  1396. for i in range(nb):
  1397. a = Reconstruction(Mactilde,i, lmonom2, [])
  1398. G2.append(a[0])
  1399. return G2
  1400.  
  1401. ###############################################
  1402. # Reduced Trop GB #
  1403. ###############################################
  1404.  
  1405. def reduced_trop_GB(F,R1,w):
  1406. r"""
  1407. Computes a reduced tropical Gröbner basis of the ideal generated by the list of polynomials in F.
  1408. The tropical term order is defined by the valuation on the base ring of R1, its ambiant
  1409. monomial order, and the tropical weight w.
  1410.  
  1411.  
  1412. INPUT:
  1413.  
  1414. - ``F`` -- a list of polynomials
  1415. - ``R1`` -- a polynomial ring
  1416. - ``w`` -- a tropical weight, i.e. an n-tuple of real numbers
  1417.  
  1418. OUTPUT:
  1419.  
  1420. A reduced tropical Gröbner basis of the ideal generated by F.
  1421.  
  1422.  
  1423. EXAMPLES::
  1424.  
  1425. sage: from sage.rings.polynomial.padics.F5Trop_doctest import *
  1426. sage: S.<x,y,z>=QQ[]
  1427. sage: R.<x1,y1,z1>=Qp(2)[]
  1428. sage: w = [0,0,0]
  1429. sage: list1 = [z^2, x^3 + 2*x^2*y + y^3, y^3 + y^2*z + y*z^2, y^2*z + z^3]
  1430. sage: G = reduced_trop_GB(list1,R,w)
  1431. sage: G
  1432. [z^2, x^3 + 2*x^2*y, y^3, y^2*z]
  1433.  
  1434. """
  1435. G = tentative_F5Trop(F,R1,w)
  1436. G = minimalization(G,R1,w)
  1437. G = reduction_trop_GB(G,R1,w)
  1438. return G
Add Comment
Please, Sign In to add comment
Public Pastes
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy.  OK, I Understand
Not a member of Pastebin yet?
Sign Up, it unlocks many cool features!