ASKSAGE: Sage Q&A Forum - RSS feedhttps://ask.sagemath.org/questions/Q&A Forum for SageenCopyright Sage, 2010. Some rights reserved under creative commons license.Fri, 09 Jul 2021 11:30:05 +0200Quick (trivial) Groebner basis but too long lift of onehttps://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/I have two (equivalent) sub-systems of polynomial equations of 10 variables and 12 equations, linked with one extra equations (so a total of 20 variable and 25 equations).
I can very quicky compute the Groebner basis of each sub-system. Then I can put them together, with the extra equation, and then I can quicky get the Groebner basis of the full system, which turns out to be trivial (total computation time: less than 3 minutes). See all the details in Appendix.
A direct computation of the Groebner basis of the full system is too long (I stopped after 1 day). The problem is that I need to compute the lift of 1 of the original full system (i.e. express 1 as a linear combination (with polynomial coefficients) of the original generators, by using the lift command `list(R.one().lift(Id))`), but it is also too long...
Now I expect that the application of a strategy involving the two sub-systems (as above) should exist to quickly get the lift of 1, but **how?**
I tried to compute this lift in two times, using the Groebner basis of the sub-systems, but it is also too long...
_____
##Appendix
First sub-system (10 variables and 12 equations)
sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
sage: A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
....: 5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
....: 25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
....: 5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
....: 5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
....: 5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
....: 5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
....: 25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
....: 5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
....: 25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
....: 5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
....: 5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
sage: Id=R.ideal(A1)
sage: %time G1=Id.groebner_basis()
CPU times: user 1.33 s, sys: 0 ns, total: 1.33 s
Wall time: 1.33 s
sage: C1=[g for g in G1]
Second (equivalent) sub-system:
sage: R.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
sage: A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
....: 5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
....: 25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
....: 5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
....: 5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
....: 5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
....: 5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
....: 25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
....: 5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
....: 5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
....: 25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
....: 5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
sage: Id=R.ideal(A2)
sage: %time G2=Id.groebner_basis()
CPU times: user 2.08 s, sys: 0 ns, total: 2.08 s
Wall time: 2.08 s
sage: C2=[g for g in G2]
The full system:
sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
sage: C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
sage: Id=R.ideal(C)
sage: %time G=Id.groebner_basis()
CPU times: user 2min 29s, sys: 16 ms, total: 2min 29s
Wall time: 2min 29s
sage: G
[1]Thu, 20 May 2021 12:28:20 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/Comment by Sébastien Palcoux for <p>I have two (equivalent) sub-systems of polynomial equations of 10 variables and 12 equations, linked with one extra equations (so a total of 20 variable and 25 equations).</p>
<p>I can very quicky compute the Groebner basis of each sub-system. Then I can put them together, with the extra equation, and then I can quicky get the Groebner basis of the full system, which turns out to be trivial (total computation time: less than 3 minutes). See all the details in Appendix.</p>
<p>A direct computation of the Groebner basis of the full system is too long (I stopped after 1 day). The problem is that I need to compute the lift of 1 of the original full system (i.e. express 1 as a linear combination (with polynomial coefficients) of the original generators, by using the lift command <code>list(R.one().lift(Id))</code>), but it is also too long... </p>
<p>Now I expect that the application of a strategy involving the two sub-systems (as above) should exist to quickly get the lift of 1, but <strong>how?</strong> </p>
<p>I tried to compute this lift in two times, using the Groebner basis of the sub-systems, but it is also too long...</p>
<hr>
<h2>Appendix</h2>
<p>First sub-system (10 variables and 12 equations)</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
sage: A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
....: 5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
....: 25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
....: 5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
....: 5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
....: 5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
....: 5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
....: 25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
....: 5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
....: 25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
....: 5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
....: 5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
sage: Id=R.ideal(A1)
sage: %time G1=Id.groebner_basis()
CPU times: user 1.33 s, sys: 0 ns, total: 1.33 s
Wall time: 1.33 s
sage: C1=[g for g in G1]
</code></pre>
<p>Second (equivalent) sub-system: </p>
<pre><code>sage: R.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
sage: A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
....: 5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
....: 25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
....: 5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
....: 5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
....: 5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
....: 5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
....: 25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
....: 5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
....: 5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
....: 25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
....: 5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
sage: Id=R.ideal(A2)
sage: %time G2=Id.groebner_basis()
CPU times: user 2.08 s, sys: 0 ns, total: 2.08 s
Wall time: 2.08 s
sage: C2=[g for g in G2]
</code></pre>
<p>The full system:</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
sage: C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
sage: Id=R.ideal(C)
sage: %time G=Id.groebner_basis()
CPU times: user 2min 29s, sys: 16 ms, total: 2min 29s
Wall time: 2min 29s
sage: G
[1]
</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57279#post-id-57279@rburing: another factor is: 321759008230308776494348949180759551. Asymptotically, the number of sqrt(n)-smooth numbers < x is known to be (1-log(2))*x + O(x/log(x)), see https://oeis.org/A063539. Now log(2)≈0.693, so there is about a 70% chance that this 785-digit number has a prime factor of more than 393 digits...Wed, 26 May 2021 08:46:53 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57279#post-id-57279Comment by rburing for <p>I have two (equivalent) sub-systems of polynomial equations of 10 variables and 12 equations, linked with one extra equations (so a total of 20 variable and 25 equations).</p>
<p>I can very quicky compute the Groebner basis of each sub-system. Then I can put them together, with the extra equation, and then I can quicky get the Groebner basis of the full system, which turns out to be trivial (total computation time: less than 3 minutes). See all the details in Appendix.</p>
<p>A direct computation of the Groebner basis of the full system is too long (I stopped after 1 day). The problem is that I need to compute the lift of 1 of the original full system (i.e. express 1 as a linear combination (with polynomial coefficients) of the original generators, by using the lift command <code>list(R.one().lift(Id))</code>), but it is also too long... </p>
<p>Now I expect that the application of a strategy involving the two sub-systems (as above) should exist to quickly get the lift of 1, but <strong>how?</strong> </p>
<p>I tried to compute this lift in two times, using the Groebner basis of the sub-systems, but it is also too long...</p>
<hr>
<h2>Appendix</h2>
<p>First sub-system (10 variables and 12 equations)</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
sage: A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
....: 5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
....: 25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
....: 5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
....: 5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
....: 5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
....: 5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
....: 25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
....: 5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
....: 25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
....: 5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
....: 5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
sage: Id=R.ideal(A1)
sage: %time G1=Id.groebner_basis()
CPU times: user 1.33 s, sys: 0 ns, total: 1.33 s
Wall time: 1.33 s
sage: C1=[g for g in G1]
</code></pre>
<p>Second (equivalent) sub-system: </p>
<pre><code>sage: R.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
sage: A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
....: 5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
....: 25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
....: 5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
....: 5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
....: 5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
....: 5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
....: 25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
....: 5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
....: 5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
....: 25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
....: 5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
sage: Id=R.ideal(A2)
sage: %time G2=Id.groebner_basis()
CPU times: user 2.08 s, sys: 0 ns, total: 2.08 s
Wall time: 2.08 s
sage: C2=[g for g in G2]
</code></pre>
<p>The full system:</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
sage: C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
sage: Id=R.ideal(C)
sage: %time G=Id.groebner_basis()
CPU times: user 2min 29s, sys: 16 ms, total: 2min 29s
Wall time: 2min 29s
sage: G
[1]
</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57261#post-id-57261I don't know what's the best way, but the factors so far were found using the [elliptic curve method](https://www.alpertron.com.ar/ECM.HTM); another factor is 6894196105010275529923254251027.Tue, 25 May 2021 09:20:05 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57261#post-id-57261Comment by Sébastien Palcoux for <p>I have two (equivalent) sub-systems of polynomial equations of 10 variables and 12 equations, linked with one extra equations (so a total of 20 variable and 25 equations).</p>
<p>I can very quicky compute the Groebner basis of each sub-system. Then I can put them together, with the extra equation, and then I can quicky get the Groebner basis of the full system, which turns out to be trivial (total computation time: less than 3 minutes). See all the details in Appendix.</p>
<p>A direct computation of the Groebner basis of the full system is too long (I stopped after 1 day). The problem is that I need to compute the lift of 1 of the original full system (i.e. express 1 as a linear combination (with polynomial coefficients) of the original generators, by using the lift command <code>list(R.one().lift(Id))</code>), but it is also too long... </p>
<p>Now I expect that the application of a strategy involving the two sub-systems (as above) should exist to quickly get the lift of 1, but <strong>how?</strong> </p>
<p>I tried to compute this lift in two times, using the Groebner basis of the sub-systems, but it is also too long...</p>
<hr>
<h2>Appendix</h2>
<p>First sub-system (10 variables and 12 equations)</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
sage: A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
....: 5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
....: 25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
....: 5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
....: 5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
....: 5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
....: 5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
....: 25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
....: 5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
....: 25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
....: 5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
....: 5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
sage: Id=R.ideal(A1)
sage: %time G1=Id.groebner_basis()
CPU times: user 1.33 s, sys: 0 ns, total: 1.33 s
Wall time: 1.33 s
sage: C1=[g for g in G1]
</code></pre>
<p>Second (equivalent) sub-system: </p>
<pre><code>sage: R.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
sage: A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
....: 5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
....: 25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
....: 5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
....: 5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
....: 5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
....: 5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
....: 25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
....: 5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
....: 5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
....: 25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
....: 5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
sage: Id=R.ideal(A2)
sage: %time G2=Id.groebner_basis()
CPU times: user 2.08 s, sys: 0 ns, total: 2.08 s
Wall time: 2.08 s
sage: C2=[g for g in G2]
</code></pre>
<p>The full system:</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
sage: C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
sage: Id=R.ideal(C)
sage: %time G=Id.groebner_basis()
CPU times: user 2min 29s, sys: 16 ms, total: 2min 29s
Wall time: 2min 29s
sage: G
[1]
</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57260#post-id-57260@rburing Which method did you use to start the factorization of this huge number? I tried the function 'trial_division' to get the next prime factor, but it did not finish after 20 hours...Tue, 25 May 2021 07:41:48 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57260#post-id-57260Comment by rburing for <p>I have two (equivalent) sub-systems of polynomial equations of 10 variables and 12 equations, linked with one extra equations (so a total of 20 variable and 25 equations).</p>
<p>I can very quicky compute the Groebner basis of each sub-system. Then I can put them together, with the extra equation, and then I can quicky get the Groebner basis of the full system, which turns out to be trivial (total computation time: less than 3 minutes). See all the details in Appendix.</p>
<p>A direct computation of the Groebner basis of the full system is too long (I stopped after 1 day). The problem is that I need to compute the lift of 1 of the original full system (i.e. express 1 as a linear combination (with polynomial coefficients) of the original generators, by using the lift command <code>list(R.one().lift(Id))</code>), but it is also too long... </p>
<p>Now I expect that the application of a strategy involving the two sub-systems (as above) should exist to quickly get the lift of 1, but <strong>how?</strong> </p>
<p>I tried to compute this lift in two times, using the Groebner basis of the sub-systems, but it is also too long...</p>
<hr>
<h2>Appendix</h2>
<p>First sub-system (10 variables and 12 equations)</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
sage: A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
....: 5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
....: 25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
....: 5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
....: 5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
....: 5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
....: 5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
....: 25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
....: 5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
....: 25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
....: 5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
....: 5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
sage: Id=R.ideal(A1)
sage: %time G1=Id.groebner_basis()
CPU times: user 1.33 s, sys: 0 ns, total: 1.33 s
Wall time: 1.33 s
sage: C1=[g for g in G1]
</code></pre>
<p>Second (equivalent) sub-system: </p>
<pre><code>sage: R.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
sage: A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
....: 5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
....: 25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
....: 5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
....: 5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
....: 5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
....: 5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
....: 25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
....: 5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
....: 5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
....: 25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
....: 5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
sage: Id=R.ideal(A2)
sage: %time G2=Id.groebner_basis()
CPU times: user 2.08 s, sys: 0 ns, total: 2.08 s
Wall time: 2.08 s
sage: C2=[g for g in G2]
</code></pre>
<p>The full system:</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
sage: C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
sage: Id=R.ideal(C)
sage: %time G=Id.groebner_basis()
CPU times: user 2min 29s, sys: 16 ms, total: 2min 29s
Wall time: 2min 29s
sage: G
[1]
</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57213#post-id-57213Numerator: `-43229237059585206766910192025630064873887338354808109523366236412917325235922459549533645895191959192066719660647997170756963396600182830269259739604438096890894083785119436868026134859982272657810992372049730997644377492574830090390702804227717772938689851430047285152773219568614271508771592335243291669125598087197566483446660761474960738664919286270417963055969225295503245692745022501334123551198196578676411370309040309672543462159719672046960825376520537565063878329462932875402735403178340699102220483445515755335938757848830865195449531261209878912968158467907562249567677739660984434298623976246460940883597744327986677370019034913909890091568068986627169835305929995631556902326395829146984403788734055514604409601849177378808235702317142670869617337832605048071523585577383`.Sat, 22 May 2021 09:47:14 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57213#post-id-57213Comment by rburing for <p>I have two (equivalent) sub-systems of polynomial equations of 10 variables and 12 equations, linked with one extra equations (so a total of 20 variable and 25 equations).</p>
<p>I can very quicky compute the Groebner basis of each sub-system. Then I can put them together, with the extra equation, and then I can quicky get the Groebner basis of the full system, which turns out to be trivial (total computation time: less than 3 minutes). See all the details in Appendix.</p>
<p>A direct computation of the Groebner basis of the full system is too long (I stopped after 1 day). The problem is that I need to compute the lift of 1 of the original full system (i.e. express 1 as a linear combination (with polynomial coefficients) of the original generators, by using the lift command <code>list(R.one().lift(Id))</code>), but it is also too long... </p>
<p>Now I expect that the application of a strategy involving the two sub-systems (as above) should exist to quickly get the lift of 1, but <strong>how?</strong> </p>
<p>I tried to compute this lift in two times, using the Groebner basis of the sub-systems, but it is also too long...</p>
<hr>
<h2>Appendix</h2>
<p>First sub-system (10 variables and 12 equations)</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
sage: A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
....: 5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
....: 25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
....: 5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
....: 5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
....: 5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
....: 5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
....: 25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
....: 5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
....: 25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
....: 5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
....: 5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
sage: Id=R.ideal(A1)
sage: %time G1=Id.groebner_basis()
CPU times: user 1.33 s, sys: 0 ns, total: 1.33 s
Wall time: 1.33 s
sage: C1=[g for g in G1]
</code></pre>
<p>Second (equivalent) sub-system: </p>
<pre><code>sage: R.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
sage: A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
....: 5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
....: 25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
....: 5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
....: 5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
....: 5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
....: 5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
....: 25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
....: 5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
....: 5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
....: 25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
....: 5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
sage: Id=R.ideal(A2)
sage: %time G2=Id.groebner_basis()
CPU times: user 2.08 s, sys: 0 ns, total: 2.08 s
Wall time: 2.08 s
sage: C2=[g for g in G2]
</code></pre>
<p>The full system:</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
sage: C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
sage: Id=R.ideal(C)
sage: %time G=Id.groebner_basis()
CPU times: user 2min 29s, sys: 16 ms, total: 2min 29s
Wall time: 2min 29s
sage: G
[1]
</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57212#post-id-57212The discrepancy is probably related to the fact that I took `C1` and `C2` which were calculated over $\mathbb{Q}$. Maybe a Groebner basis over $\mathbb{Z}$ (after clearing denominators of generators) would work better? I don't know much about this or about characteristic $p$. The denominator of the determinant is $2^{560} \cdot 3^{343} \cdot 5^{588} \cdot 7^{198} \cdot 19^{106} \cdot 29^{83}$.Sat, 22 May 2021 09:46:08 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57212#post-id-57212Comment by Sébastien Palcoux for <p>I have two (equivalent) sub-systems of polynomial equations of 10 variables and 12 equations, linked with one extra equations (so a total of 20 variable and 25 equations).</p>
<p>I can very quicky compute the Groebner basis of each sub-system. Then I can put them together, with the extra equation, and then I can quicky get the Groebner basis of the full system, which turns out to be trivial (total computation time: less than 3 minutes). See all the details in Appendix.</p>
<p>A direct computation of the Groebner basis of the full system is too long (I stopped after 1 day). The problem is that I need to compute the lift of 1 of the original full system (i.e. express 1 as a linear combination (with polynomial coefficients) of the original generators, by using the lift command <code>list(R.one().lift(Id))</code>), but it is also too long... </p>
<p>Now I expect that the application of a strategy involving the two sub-systems (as above) should exist to quickly get the lift of 1, but <strong>how?</strong> </p>
<p>I tried to compute this lift in two times, using the Groebner basis of the sub-systems, but it is also too long...</p>
<hr>
<h2>Appendix</h2>
<p>First sub-system (10 variables and 12 equations)</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
sage: A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
....: 5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
....: 25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
....: 5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
....: 5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
....: 5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
....: 5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
....: 25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
....: 5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
....: 25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
....: 5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
....: 5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
sage: Id=R.ideal(A1)
sage: %time G1=Id.groebner_basis()
CPU times: user 1.33 s, sys: 0 ns, total: 1.33 s
Wall time: 1.33 s
sage: C1=[g for g in G1]
</code></pre>
<p>Second (equivalent) sub-system: </p>
<pre><code>sage: R.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
sage: A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
....: 5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
....: 25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
....: 5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
....: 5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
....: 5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
....: 5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
....: 25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
....: 5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
....: 5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
....: 25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
....: 5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
sage: Id=R.ideal(A2)
sage: %time G2=Id.groebner_basis()
CPU times: user 2.08 s, sys: 0 ns, total: 2.08 s
Wall time: 2.08 s
sage: C2=[g for g in G2]
</code></pre>
<p>The full system:</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
sage: C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
sage: Id=R.ideal(C)
sage: %time G=Id.groebner_basis()
CPU times: user 2min 29s, sys: 16 ms, total: 2min 29s
Wall time: 2min 29s
sage: G
[1]
</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57206#post-id-57206@rburing Very interesting!! If I am not mistaken, the primes that I am looking for are exactly the prime divisors of the determinant of this 196x196 matrix, right?
If so, 29 and 41 should also be divisors. Did you miss them, or am I mistaken?
Can you write the numerator of the determinant (with its 785 digits)?Fri, 21 May 2021 21:47:23 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57206#post-id-57206Comment by rburing for <p>I have two (equivalent) sub-systems of polynomial equations of 10 variables and 12 equations, linked with one extra equations (so a total of 20 variable and 25 equations).</p>
<p>I can very quicky compute the Groebner basis of each sub-system. Then I can put them together, with the extra equation, and then I can quicky get the Groebner basis of the full system, which turns out to be trivial (total computation time: less than 3 minutes). See all the details in Appendix.</p>
<p>A direct computation of the Groebner basis of the full system is too long (I stopped after 1 day). The problem is that I need to compute the lift of 1 of the original full system (i.e. express 1 as a linear combination (with polynomial coefficients) of the original generators, by using the lift command <code>list(R.one().lift(Id))</code>), but it is also too long... </p>
<p>Now I expect that the application of a strategy involving the two sub-systems (as above) should exist to quickly get the lift of 1, but <strong>how?</strong> </p>
<p>I tried to compute this lift in two times, using the Groebner basis of the sub-systems, but it is also too long...</p>
<hr>
<h2>Appendix</h2>
<p>First sub-system (10 variables and 12 equations)</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
sage: A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
....: 5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
....: 25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
....: 5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
....: 5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
....: 5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
....: 5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
....: 25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
....: 5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
....: 25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
....: 5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
....: 5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
sage: Id=R.ideal(A1)
sage: %time G1=Id.groebner_basis()
CPU times: user 1.33 s, sys: 0 ns, total: 1.33 s
Wall time: 1.33 s
sage: C1=[g for g in G1]
</code></pre>
<p>Second (equivalent) sub-system: </p>
<pre><code>sage: R.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
sage: A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
....: 5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
....: 25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
....: 5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
....: 5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
....: 5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
....: 5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
....: 25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
....: 5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
....: 5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
....: 25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
....: 5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
sage: Id=R.ideal(A2)
sage: %time G2=Id.groebner_basis()
CPU times: user 2.08 s, sys: 0 ns, total: 2.08 s
Wall time: 2.08 s
sage: C2=[g for g in G2]
</code></pre>
<p>The full system:</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
sage: C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
sage: Id=R.ideal(C)
sage: %time G=Id.groebner_basis()
CPU times: user 2min 29s, sys: 16 ms, total: 2min 29s
Wall time: 2min 29s
sage: G
[1]
</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57198#post-id-57198Interesting problem. The numerator of the determinant of the aforementioned matrix has 785 digits and its factorization is 11 × 467 × 1423 × 242057 × 1477913 × 2673859 × 9528289 × 3869785541 × 715114922077295989 × 75154031037086343317 × 652748578372455171701 × (to be determined).Fri, 21 May 2021 12:37:33 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57198#post-id-57198Comment by Sébastien Palcoux for <p>I have two (equivalent) sub-systems of polynomial equations of 10 variables and 12 equations, linked with one extra equations (so a total of 20 variable and 25 equations).</p>
<p>I can very quicky compute the Groebner basis of each sub-system. Then I can put them together, with the extra equation, and then I can quicky get the Groebner basis of the full system, which turns out to be trivial (total computation time: less than 3 minutes). See all the details in Appendix.</p>
<p>A direct computation of the Groebner basis of the full system is too long (I stopped after 1 day). The problem is that I need to compute the lift of 1 of the original full system (i.e. express 1 as a linear combination (with polynomial coefficients) of the original generators, by using the lift command <code>list(R.one().lift(Id))</code>), but it is also too long... </p>
<p>Now I expect that the application of a strategy involving the two sub-systems (as above) should exist to quickly get the lift of 1, but <strong>how?</strong> </p>
<p>I tried to compute this lift in two times, using the Groebner basis of the sub-systems, but it is also too long...</p>
<hr>
<h2>Appendix</h2>
<p>First sub-system (10 variables and 12 equations)</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
sage: A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
....: 5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
....: 25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
....: 5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
....: 5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
....: 5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
....: 5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
....: 25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
....: 5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
....: 25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
....: 5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
....: 5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
sage: Id=R.ideal(A1)
sage: %time G1=Id.groebner_basis()
CPU times: user 1.33 s, sys: 0 ns, total: 1.33 s
Wall time: 1.33 s
sage: C1=[g for g in G1]
</code></pre>
<p>Second (equivalent) sub-system: </p>
<pre><code>sage: R.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
sage: A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
....: 5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
....: 25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
....: 5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
....: 5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
....: 5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
....: 5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
....: 25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
....: 5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
....: 5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
....: 25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
....: 5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
sage: Id=R.ideal(A2)
sage: %time G2=Id.groebner_basis()
CPU times: user 2.08 s, sys: 0 ns, total: 2.08 s
Wall time: 2.08 s
sage: C2=[g for g in G2]
</code></pre>
<p>The full system:</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
sage: C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
sage: Id=R.ideal(C)
sage: %time G=Id.groebner_basis()
CPU times: user 2min 29s, sys: 16 ms, total: 2min 29s
Wall time: 2min 29s
sage: G
[1]
</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57192#post-id-57192@rburing What you wrote in your first comment was done firstly, but that requires to use numerical solutions... The computation with Groebner bases only is cleaner (in my opinion) and still very quick in this case. Anyway, my problem is that I need the lift of 1, because the Groebner basis must be also trivial in char p>0, except for finitely many p, contained in the set of prime divisors of the denominators in the lift of 1, but I need to classify all p for which the Groebner basis is non-trivial: for p<1000000, there are exactly six one: 11, 29, 41, 467, 1423, 242057. With a lift of 1, I can easily classify all such p.Thu, 20 May 2021 20:53:56 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57192#post-id-57192Comment by rburing for <p>I have two (equivalent) sub-systems of polynomial equations of 10 variables and 12 equations, linked with one extra equations (so a total of 20 variable and 25 equations).</p>
<p>I can very quicky compute the Groebner basis of each sub-system. Then I can put them together, with the extra equation, and then I can quicky get the Groebner basis of the full system, which turns out to be trivial (total computation time: less than 3 minutes). See all the details in Appendix.</p>
<p>A direct computation of the Groebner basis of the full system is too long (I stopped after 1 day). The problem is that I need to compute the lift of 1 of the original full system (i.e. express 1 as a linear combination (with polynomial coefficients) of the original generators, by using the lift command <code>list(R.one().lift(Id))</code>), but it is also too long... </p>
<p>Now I expect that the application of a strategy involving the two sub-systems (as above) should exist to quickly get the lift of 1, but <strong>how?</strong> </p>
<p>I tried to compute this lift in two times, using the Groebner basis of the sub-systems, but it is also too long...</p>
<hr>
<h2>Appendix</h2>
<p>First sub-system (10 variables and 12 equations)</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
sage: A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
....: 5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
....: 25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
....: 5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
....: 5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
....: 5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
....: 5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
....: 25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
....: 5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
....: 25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
....: 5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
....: 5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
sage: Id=R.ideal(A1)
sage: %time G1=Id.groebner_basis()
CPU times: user 1.33 s, sys: 0 ns, total: 1.33 s
Wall time: 1.33 s
sage: C1=[g for g in G1]
</code></pre>
<p>Second (equivalent) sub-system: </p>
<pre><code>sage: R.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
sage: A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
....: 5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
....: 25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
....: 5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
....: 5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
....: 5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
....: 5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
....: 25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
....: 5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
....: 5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
....: 25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
....: 5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
sage: Id=R.ideal(A2)
sage: %time G2=Id.groebner_basis()
CPU times: user 2.08 s, sys: 0 ns, total: 2.08 s
Wall time: 2.08 s
sage: C2=[g for g in G2]
</code></pre>
<p>The full system:</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
sage: C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
sage: Id=R.ideal(C)
sage: %time G=Id.groebner_basis()
CPU times: user 2min 29s, sys: 16 ms, total: 2min 29s
Wall time: 2min 29s
sage: G
[1]
</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57191#post-id-57191A simpler variant of that idea (using linear algebra over $\mathbb{Q}$ instead of arithmetic over $\bar{\mathbb{Q}}$) is to use that `C1 + C2` is a Groebner basis of the ideal $I$ generated by `C1 + C2`; hence there is a normal basis (basis of the quotient $R/I$ as a $\mathbb{Q}$-vector space) consisting of 196 monomials, and multiplication by `f = 5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125` can be written as a 196x196 matrix over $\mathbb{Q}$, and the values of `f` on $V(I)$ are the eigenvalues of this matrix, so it suffices to show that the matrix does not have 0 as an eigenvalue, or that it has full rank, or that its determinant is nonzero (which is easy).Thu, 20 May 2021 20:27:03 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57191#post-id-57191Comment by rburing for <p>I have two (equivalent) sub-systems of polynomial equations of 10 variables and 12 equations, linked with one extra equations (so a total of 20 variable and 25 equations).</p>
<p>I can very quicky compute the Groebner basis of each sub-system. Then I can put them together, with the extra equation, and then I can quicky get the Groebner basis of the full system, which turns out to be trivial (total computation time: less than 3 minutes). See all the details in Appendix.</p>
<p>A direct computation of the Groebner basis of the full system is too long (I stopped after 1 day). The problem is that I need to compute the lift of 1 of the original full system (i.e. express 1 as a linear combination (with polynomial coefficients) of the original generators, by using the lift command <code>list(R.one().lift(Id))</code>), but it is also too long... </p>
<p>Now I expect that the application of a strategy involving the two sub-systems (as above) should exist to quickly get the lift of 1, but <strong>how?</strong> </p>
<p>I tried to compute this lift in two times, using the Groebner basis of the sub-systems, but it is also too long...</p>
<hr>
<h2>Appendix</h2>
<p>First sub-system (10 variables and 12 equations)</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
sage: A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
....: 5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
....: 25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
....: 5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
....: 5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
....: 5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
....: 5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
....: 25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
....: 5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
....: 25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
....: 5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
....: 5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
sage: Id=R.ideal(A1)
sage: %time G1=Id.groebner_basis()
CPU times: user 1.33 s, sys: 0 ns, total: 1.33 s
Wall time: 1.33 s
sage: C1=[g for g in G1]
</code></pre>
<p>Second (equivalent) sub-system: </p>
<pre><code>sage: R.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
sage: A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
....: 5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
....: 25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
....: 5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
....: 5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
....: 5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
....: 5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
....: 25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
....: 5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
....: 5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
....: 25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
....: 5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
sage: Id=R.ideal(A2)
sage: %time G2=Id.groebner_basis()
CPU times: user 2.08 s, sys: 0 ns, total: 2.08 s
Wall time: 2.08 s
sage: C2=[g for g in G2]
</code></pre>
<p>The full system:</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
sage: C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
sage: Id=R.ideal(C)
sage: %time G=Id.groebner_basis()
CPU times: user 2min 29s, sys: 16 ms, total: 2min 29s
Wall time: 2min 29s
sage: G
[1]
</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57183#post-id-57183It's not what you asked, and probably suboptimal, but: since the two individual ideals are zero-dimensional, and the respective varieties each contain 14 points over $\bar{\mathbb{Q}}$, you could also loop through all 196 pairs of points and substitute them into the extra polynomial, and check that you never get zero.Thu, 20 May 2021 16:00:32 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57183#post-id-57183Answer by rburing for <p>I have two (equivalent) sub-systems of polynomial equations of 10 variables and 12 equations, linked with one extra equations (so a total of 20 variable and 25 equations).</p>
<p>I can very quicky compute the Groebner basis of each sub-system. Then I can put them together, with the extra equation, and then I can quicky get the Groebner basis of the full system, which turns out to be trivial (total computation time: less than 3 minutes). See all the details in Appendix.</p>
<p>A direct computation of the Groebner basis of the full system is too long (I stopped after 1 day). The problem is that I need to compute the lift of 1 of the original full system (i.e. express 1 as a linear combination (with polynomial coefficients) of the original generators, by using the lift command <code>list(R.one().lift(Id))</code>), but it is also too long... </p>
<p>Now I expect that the application of a strategy involving the two sub-systems (as above) should exist to quickly get the lift of 1, but <strong>how?</strong> </p>
<p>I tried to compute this lift in two times, using the Groebner basis of the sub-systems, but it is also too long...</p>
<hr>
<h2>Appendix</h2>
<p>First sub-system (10 variables and 12 equations)</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
sage: A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
....: 5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
....: 25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
....: 5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
....: 5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
....: 5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
....: 5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
....: 25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
....: 5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
....: 25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
....: 5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
....: 5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
sage: Id=R.ideal(A1)
sage: %time G1=Id.groebner_basis()
CPU times: user 1.33 s, sys: 0 ns, total: 1.33 s
Wall time: 1.33 s
sage: C1=[g for g in G1]
</code></pre>
<p>Second (equivalent) sub-system: </p>
<pre><code>sage: R.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
sage: A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
....: 5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
....: 25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
....: 5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
....: 5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
....: 5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
....: 5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
....: 25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
....: 5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
....: 5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
....: 25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
....: 5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
sage: Id=R.ideal(A2)
sage: %time G2=Id.groebner_basis()
CPU times: user 2.08 s, sys: 0 ns, total: 2.08 s
Wall time: 2.08 s
sage: C2=[g for g in G2]
</code></pre>
<p>The full system:</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
sage: C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
sage: Id=R.ideal(C)
sage: %time G=Id.groebner_basis()
CPU times: user 2min 29s, sys: 16 ms, total: 2min 29s
Wall time: 2min 29s
sage: G
[1]
</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?answer=57244#post-id-57244In this rather special case we can find a lift of $1$ quickly and efficiently using linear algebra.
Let $R$ be a polynomial ring over $\mathbb{C}$, $I \subset R$ a zero-dimensional radical ideal and $f \in R$. To show $ \langle 1 \rangle = I + \langle f \rangle$ it suffices (by Hilbert's Nullstellensatz) to show $\varnothing = V(I + \langle f \rangle) = V(I) \cap V(f)$, i.e. that $f$ does not vanish on $V(I)$. Since $I$ is radical, the ideal of polynomials in $R$ vanishing on $V(I)$ is exactly $I$ (by the Nullstellensatz). By definition of zero-dimensional, the quotient $R/I$ is a finite-dimensional $\mathbb{C}$-vector space.
* Consider the $\mathbb{C}$-linear "multiplication by $f$" map $M_f: g \mapsto f \cdot g$ on $R/I$. Let $g$ be an eigenvector of $M_f$ with eigenvalue $\lambda$, so $f \cdot g = \lambda \cdot g$ or $(f-\lambda) g = 0$ holds in $R/I$, i.e. $(f-\lambda)g \in I$. The fact that $g \neq 0$ in $R/I$ means $g \not\in I$, hence there is a point $p \in V(I)$ with $g(p) \neq 0$. Now from $(f-\lambda)g \in I$ it follows that $\lambda = f(p)$ is the value of $f$ at the point $p$.
* Conversely, if $P$ is the minimum polynomial of $M_f$ then $P(M_f) = 0$ implies $P(f) = 0$ in $R/I$ (by evaluating at $1 \in R/I$), hence $P(f) \in I$ and it follows that $0 = P(f)(p) = P(f(p))$ for all $p \in V(I)$, i.e. all values $f(p)$ for $p \in V(I)$ are zeros of $P$, and hence they are eigenvalues of $M_f$.
So the values of $f$ on $V(I)$ are exactly the eigenvalues of $M_f$, and it suffices to show that $0$ is not one of them.
Given a Groebner basis $G$ of $I$, the multivariate polynomial division algorithm implies that the set of all monomials not divisible by any leading terms of $G$ is a basis of $R/I$ as a $\mathbb{C}$-vector space (called the normal basis). The matrix of $M_f$ with respect to this basis is easily computed:
def multiplication_matrix(f, I):
G = I.groebner_basis()
f = f.reduce(G)
NB = I.normal_basis()
M = Matrix(f.base_ring(), len(NB))
for (i, m1) in enumerate(NB):
g = (f*m1).reduce(G)
M.set_column(i, [g.monomial_coefficient(m2) for m2 in NB])
return M
We apply it to the problem at hand:
sage: f = 5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125
sage: I = R.ideal(C1 + C2)
sage: M_f = multiplication_matrix(f, I)
sage: M_f.det() != 0
True
Hence $f$ does not attain the value $0$ on $V(I)$, so $\varnothing = V(I + \langle f \rangle)$ and $\langle 1 \rangle = I + \langle f \rangle$.
Since $M_f$ is invertible, we can find a polynomial representative of $f^{-1} + I$:
sage: NB = list(I.normal_basis())
sage: len(NB), NB.index(1)
(196, 195)
sage: finv_coeffs = M_f \ vector([0]*195 + [1])
sage: finv = sum(c*m for (c,m) in zip(finv_coeffs, NB))
sage: (f*finv).reduce(I)
1
Then we have $1 = f^{-1} \cdot f$ in $R/I$, i.e. $1 = f^{-1} \cdot f + g$ for some $g \in I$, and we can express $g$ in terms of generators by applying the division algorithm to $1 - f^{-1} \cdot f$ (or finding a lift of it to $I$):
sage: lift_coeffs = [finv] + list((1 - finv*f).lift(I))
sage: sum(c*m for c,m in zip(lift_coeffs, [f] + I.gens()))
1Mon, 24 May 2021 01:25:28 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?answer=57244#post-id-57244Comment by Sébastien Palcoux for <p>In this rather special case we can find a lift of $1$ quickly and efficiently using linear algebra.</p>
<p>Let $R$ be a polynomial ring over $\mathbb{C}$, $I \subset R$ a zero-dimensional radical ideal and $f \in R$. To show $ \langle 1 \rangle = I + \langle f \rangle$ it suffices (by Hilbert's Nullstellensatz) to show $\varnothing = V(I + \langle f \rangle) = V(I) \cap V(f)$, i.e. that $f$ does not vanish on $V(I)$. Since $I$ is radical, the ideal of polynomials in $R$ vanishing on $V(I)$ is exactly $I$ (by the Nullstellensatz). By definition of zero-dimensional, the quotient $R/I$ is a finite-dimensional $\mathbb{C}$-vector space.</p>
<ul>
<li><p>Consider the $\mathbb{C}$-linear "multiplication by $f$" map $M_f: g \mapsto f \cdot g$ on $R/I$. Let $g$ be an eigenvector of $M_f$ with eigenvalue $\lambda$, so $f \cdot g = \lambda \cdot g$ or $(f-\lambda) g = 0$ holds in $R/I$, i.e. $(f-\lambda)g \in I$. The fact that $g \neq 0$ in $R/I$ means $g \not\in I$, hence there is a point $p \in V(I)$ with $g(p) \neq 0$. Now from $(f-\lambda)g \in I$ it follows that $\lambda = f(p)$ is the value of $f$ at the point $p$.</p></li>
<li><p>Conversely, if $P$ is the minimum polynomial of $M_f$ then $P(M_f) = 0$ implies $P(f) = 0$ in $R/I$ (by evaluating at $1 \in R/I$), hence $P(f) \in I$ and it follows that $0 = P(f)(p) = P(f(p))$ for all $p \in V(I)$, i.e. all values $f(p)$ for $p \in V(I)$ are zeros of $P$, and hence they are eigenvalues of $M_f$.</p></li>
</ul>
<p>So the values of $f$ on $V(I)$ are exactly the eigenvalues of $M_f$, and it suffices to show that $0$ is not one of them.</p>
<p>Given a Groebner basis $G$ of $I$, the multivariate polynomial division algorithm implies that the set of all monomials not divisible by any leading terms of $G$ is a basis of $R/I$ as a $\mathbb{C}$-vector space (called the normal basis). The matrix of $M_f$ with respect to this basis is easily computed:</p>
<pre><code>def multiplication_matrix(f, I):
G = I.groebner_basis()
f = f.reduce(G)
NB = I.normal_basis()
M = Matrix(f.base_ring(), len(NB))
for (i, m1) in enumerate(NB):
g = (f*m1).reduce(G)
M.set_column(i, [g.monomial_coefficient(m2) for m2 in NB])
return M
</code></pre>
<p>We apply it to the problem at hand:</p>
<pre><code>sage: f = 5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125
sage: I = R.ideal(C1 + C2)
sage: M_f = multiplication_matrix(f, I)
sage: M_f.det() != 0
True
</code></pre>
<p>Hence $f$ does not attain the value $0$ on $V(I)$, so $\varnothing = V(I + \langle f \rangle)$ and $\langle 1 \rangle = I + \langle f \rangle$.</p>
<p>Since $M_f$ is invertible, we can find a polynomial representative of $f^{-1} + I$:</p>
<pre><code>sage: NB = list(I.normal_basis())
sage: len(NB), NB.index(1)
(196, 195)
sage: finv_coeffs = M_f \ vector([0]*195 + [1])
sage: finv = sum(c*m for (c,m) in zip(finv_coeffs, NB))
sage: (f*finv).reduce(I)
1
</code></pre>
<p>Then we have $1 = f^{-1} \cdot f$ in $R/I$, i.e. $1 = f^{-1} \cdot f + g$ for some $g \in I$, and we can express $g$ in terms of generators by applying the division algorithm to $1 - f^{-1} \cdot f$ (or finding a lift of it to $I$):</p>
<pre><code>sage: lift_coeffs = [finv] + list((1 - finv*f).lift(I))
sage: sum(c*m for c,m in zip(lift_coeffs, [f] + I.gens()))
1
</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57936#post-id-57936@rburing see [this new post](https://ask.sagemath.org/question/57934/polynomial-system-without-solution-in-char-0-classification-of-char-p-with-solution/) motivating the problem in my previous comment.Fri, 09 Jul 2021 11:30:05 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57936#post-id-57936Comment by Sébastien Palcoux for <p>In this rather special case we can find a lift of $1$ quickly and efficiently using linear algebra.</p>
<p>Let $R$ be a polynomial ring over $\mathbb{C}$, $I \subset R$ a zero-dimensional radical ideal and $f \in R$. To show $ \langle 1 \rangle = I + \langle f \rangle$ it suffices (by Hilbert's Nullstellensatz) to show $\varnothing = V(I + \langle f \rangle) = V(I) \cap V(f)$, i.e. that $f$ does not vanish on $V(I)$. Since $I$ is radical, the ideal of polynomials in $R$ vanishing on $V(I)$ is exactly $I$ (by the Nullstellensatz). By definition of zero-dimensional, the quotient $R/I$ is a finite-dimensional $\mathbb{C}$-vector space.</p>
<ul>
<li><p>Consider the $\mathbb{C}$-linear "multiplication by $f$" map $M_f: g \mapsto f \cdot g$ on $R/I$. Let $g$ be an eigenvector of $M_f$ with eigenvalue $\lambda$, so $f \cdot g = \lambda \cdot g$ or $(f-\lambda) g = 0$ holds in $R/I$, i.e. $(f-\lambda)g \in I$. The fact that $g \neq 0$ in $R/I$ means $g \not\in I$, hence there is a point $p \in V(I)$ with $g(p) \neq 0$. Now from $(f-\lambda)g \in I$ it follows that $\lambda = f(p)$ is the value of $f$ at the point $p$.</p></li>
<li><p>Conversely, if $P$ is the minimum polynomial of $M_f$ then $P(M_f) = 0$ implies $P(f) = 0$ in $R/I$ (by evaluating at $1 \in R/I$), hence $P(f) \in I$ and it follows that $0 = P(f)(p) = P(f(p))$ for all $p \in V(I)$, i.e. all values $f(p)$ for $p \in V(I)$ are zeros of $P$, and hence they are eigenvalues of $M_f$.</p></li>
</ul>
<p>So the values of $f$ on $V(I)$ are exactly the eigenvalues of $M_f$, and it suffices to show that $0$ is not one of them.</p>
<p>Given a Groebner basis $G$ of $I$, the multivariate polynomial division algorithm implies that the set of all monomials not divisible by any leading terms of $G$ is a basis of $R/I$ as a $\mathbb{C}$-vector space (called the normal basis). The matrix of $M_f$ with respect to this basis is easily computed:</p>
<pre><code>def multiplication_matrix(f, I):
G = I.groebner_basis()
f = f.reduce(G)
NB = I.normal_basis()
M = Matrix(f.base_ring(), len(NB))
for (i, m1) in enumerate(NB):
g = (f*m1).reduce(G)
M.set_column(i, [g.monomial_coefficient(m2) for m2 in NB])
return M
</code></pre>
<p>We apply it to the problem at hand:</p>
<pre><code>sage: f = 5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125
sage: I = R.ideal(C1 + C2)
sage: M_f = multiplication_matrix(f, I)
sage: M_f.det() != 0
True
</code></pre>
<p>Hence $f$ does not attain the value $0$ on $V(I)$, so $\varnothing = V(I + \langle f \rangle)$ and $\langle 1 \rangle = I + \langle f \rangle$.</p>
<p>Since $M_f$ is invertible, we can find a polynomial representative of $f^{-1} + I$:</p>
<pre><code>sage: NB = list(I.normal_basis())
sage: len(NB), NB.index(1)
(196, 195)
sage: finv_coeffs = M_f \ vector([0]*195 + [1])
sage: finv = sum(c*m for (c,m) in zip(finv_coeffs, NB))
sage: (f*finv).reduce(I)
1
</code></pre>
<p>Then we have $1 = f^{-1} \cdot f$ in $R/I$, i.e. $1 = f^{-1} \cdot f + g$ for some $g \in I$, and we can express $g$ in terms of generators by applying the division algorithm to $1 - f^{-1} \cdot f$ (or finding a lift of it to $I$):</p>
<pre><code>sage: lift_coeffs = [finv] + list((1 - finv*f).lift(I))
sage: sum(c*m for c,m in zip(lift_coeffs, [f] + I.gens()))
1
</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57851#post-id-57851@rburing Can your method by adapted to lift $C_1$ with respect to $A_1$?Fri, 02 Jul 2021 15:20:30 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57851#post-id-57851Comment by mwageringel for <p>In this rather special case we can find a lift of $1$ quickly and efficiently using linear algebra.</p>
<p>Let $R$ be a polynomial ring over $\mathbb{C}$, $I \subset R$ a zero-dimensional radical ideal and $f \in R$. To show $ \langle 1 \rangle = I + \langle f \rangle$ it suffices (by Hilbert's Nullstellensatz) to show $\varnothing = V(I + \langle f \rangle) = V(I) \cap V(f)$, i.e. that $f$ does not vanish on $V(I)$. Since $I$ is radical, the ideal of polynomials in $R$ vanishing on $V(I)$ is exactly $I$ (by the Nullstellensatz). By definition of zero-dimensional, the quotient $R/I$ is a finite-dimensional $\mathbb{C}$-vector space.</p>
<ul>
<li><p>Consider the $\mathbb{C}$-linear "multiplication by $f$" map $M_f: g \mapsto f \cdot g$ on $R/I$. Let $g$ be an eigenvector of $M_f$ with eigenvalue $\lambda$, so $f \cdot g = \lambda \cdot g$ or $(f-\lambda) g = 0$ holds in $R/I$, i.e. $(f-\lambda)g \in I$. The fact that $g \neq 0$ in $R/I$ means $g \not\in I$, hence there is a point $p \in V(I)$ with $g(p) \neq 0$. Now from $(f-\lambda)g \in I$ it follows that $\lambda = f(p)$ is the value of $f$ at the point $p$.</p></li>
<li><p>Conversely, if $P$ is the minimum polynomial of $M_f$ then $P(M_f) = 0$ implies $P(f) = 0$ in $R/I$ (by evaluating at $1 \in R/I$), hence $P(f) \in I$ and it follows that $0 = P(f)(p) = P(f(p))$ for all $p \in V(I)$, i.e. all values $f(p)$ for $p \in V(I)$ are zeros of $P$, and hence they are eigenvalues of $M_f$.</p></li>
</ul>
<p>So the values of $f$ on $V(I)$ are exactly the eigenvalues of $M_f$, and it suffices to show that $0$ is not one of them.</p>
<p>Given a Groebner basis $G$ of $I$, the multivariate polynomial division algorithm implies that the set of all monomials not divisible by any leading terms of $G$ is a basis of $R/I$ as a $\mathbb{C}$-vector space (called the normal basis). The matrix of $M_f$ with respect to this basis is easily computed:</p>
<pre><code>def multiplication_matrix(f, I):
G = I.groebner_basis()
f = f.reduce(G)
NB = I.normal_basis()
M = Matrix(f.base_ring(), len(NB))
for (i, m1) in enumerate(NB):
g = (f*m1).reduce(G)
M.set_column(i, [g.monomial_coefficient(m2) for m2 in NB])
return M
</code></pre>
<p>We apply it to the problem at hand:</p>
<pre><code>sage: f = 5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125
sage: I = R.ideal(C1 + C2)
sage: M_f = multiplication_matrix(f, I)
sage: M_f.det() != 0
True
</code></pre>
<p>Hence $f$ does not attain the value $0$ on $V(I)$, so $\varnothing = V(I + \langle f \rangle)$ and $\langle 1 \rangle = I + \langle f \rangle$.</p>
<p>Since $M_f$ is invertible, we can find a polynomial representative of $f^{-1} + I$:</p>
<pre><code>sage: NB = list(I.normal_basis())
sage: len(NB), NB.index(1)
(196, 195)
sage: finv_coeffs = M_f \ vector([0]*195 + [1])
sage: finv = sum(c*m for (c,m) in zip(finv_coeffs, NB))
sage: (f*finv).reduce(I)
1
</code></pre>
<p>Then we have $1 = f^{-1} \cdot f$ in $R/I$, i.e. $1 = f^{-1} \cdot f + g$ for some $g \in I$, and we can express $g$ in terms of generators by applying the division algorithm to $1 - f^{-1} \cdot f$ (or finding a lift of it to $I$):</p>
<pre><code>sage: lift_coeffs = [finv] + list((1 - finv*f).lift(I))
sage: sum(c*m for c,m in zip(lift_coeffs, [f] + I.gens()))
1
</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57252#post-id-57252@rburing This is a very nice answer.Mon, 24 May 2021 11:06:41 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57252#post-id-57252Comment by Sébastien Palcoux for <p>In this rather special case we can find a lift of $1$ quickly and efficiently using linear algebra.</p>
<p>Let $R$ be a polynomial ring over $\mathbb{C}$, $I \subset R$ a zero-dimensional radical ideal and $f \in R$. To show $ \langle 1 \rangle = I + \langle f \rangle$ it suffices (by Hilbert's Nullstellensatz) to show $\varnothing = V(I + \langle f \rangle) = V(I) \cap V(f)$, i.e. that $f$ does not vanish on $V(I)$. Since $I$ is radical, the ideal of polynomials in $R$ vanishing on $V(I)$ is exactly $I$ (by the Nullstellensatz). By definition of zero-dimensional, the quotient $R/I$ is a finite-dimensional $\mathbb{C}$-vector space.</p>
<ul>
<li><p>Consider the $\mathbb{C}$-linear "multiplication by $f$" map $M_f: g \mapsto f \cdot g$ on $R/I$. Let $g$ be an eigenvector of $M_f$ with eigenvalue $\lambda$, so $f \cdot g = \lambda \cdot g$ or $(f-\lambda) g = 0$ holds in $R/I$, i.e. $(f-\lambda)g \in I$. The fact that $g \neq 0$ in $R/I$ means $g \not\in I$, hence there is a point $p \in V(I)$ with $g(p) \neq 0$. Now from $(f-\lambda)g \in I$ it follows that $\lambda = f(p)$ is the value of $f$ at the point $p$.</p></li>
<li><p>Conversely, if $P$ is the minimum polynomial of $M_f$ then $P(M_f) = 0$ implies $P(f) = 0$ in $R/I$ (by evaluating at $1 \in R/I$), hence $P(f) \in I$ and it follows that $0 = P(f)(p) = P(f(p))$ for all $p \in V(I)$, i.e. all values $f(p)$ for $p \in V(I)$ are zeros of $P$, and hence they are eigenvalues of $M_f$.</p></li>
</ul>
<p>So the values of $f$ on $V(I)$ are exactly the eigenvalues of $M_f$, and it suffices to show that $0$ is not one of them.</p>
<p>Given a Groebner basis $G$ of $I$, the multivariate polynomial division algorithm implies that the set of all monomials not divisible by any leading terms of $G$ is a basis of $R/I$ as a $\mathbb{C}$-vector space (called the normal basis). The matrix of $M_f$ with respect to this basis is easily computed:</p>
<pre><code>def multiplication_matrix(f, I):
G = I.groebner_basis()
f = f.reduce(G)
NB = I.normal_basis()
M = Matrix(f.base_ring(), len(NB))
for (i, m1) in enumerate(NB):
g = (f*m1).reduce(G)
M.set_column(i, [g.monomial_coefficient(m2) for m2 in NB])
return M
</code></pre>
<p>We apply it to the problem at hand:</p>
<pre><code>sage: f = 5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125
sage: I = R.ideal(C1 + C2)
sage: M_f = multiplication_matrix(f, I)
sage: M_f.det() != 0
True
</code></pre>
<p>Hence $f$ does not attain the value $0$ on $V(I)$, so $\varnothing = V(I + \langle f \rangle)$ and $\langle 1 \rangle = I + \langle f \rangle$.</p>
<p>Since $M_f$ is invertible, we can find a polynomial representative of $f^{-1} + I$:</p>
<pre><code>sage: NB = list(I.normal_basis())
sage: len(NB), NB.index(1)
(196, 195)
sage: finv_coeffs = M_f \ vector([0]*195 + [1])
sage: finv = sum(c*m for (c,m) in zip(finv_coeffs, NB))
sage: (f*finv).reduce(I)
1
</code></pre>
<p>Then we have $1 = f^{-1} \cdot f$ in $R/I$, i.e. $1 = f^{-1} \cdot f + g$ for some $g \in I$, and we can express $g$ in terms of generators by applying the division algorithm to $1 - f^{-1} \cdot f$ (or finding a lift of it to $I$):</p>
<pre><code>sage: lift_coeffs = [finv] + list((1 - finv*f).lift(I))
sage: sum(c*m for c,m in zip(lift_coeffs, [f] + I.gens()))
1
</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57250#post-id-57250What your answer give is the lift of $1$ with respect to $C_1 + C_2 + $ extra. What is asked is the lift of $1$ with respect to $A_1 + A_2 + $ extra (with a bit of luck, it could be easier...).Mon, 24 May 2021 10:47:25 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57250#post-id-57250Comment by Sébastien Palcoux for <p>In this rather special case we can find a lift of $1$ quickly and efficiently using linear algebra.</p>
<p>Let $R$ be a polynomial ring over $\mathbb{C}$, $I \subset R$ a zero-dimensional radical ideal and $f \in R$. To show $ \langle 1 \rangle = I + \langle f \rangle$ it suffices (by Hilbert's Nullstellensatz) to show $\varnothing = V(I + \langle f \rangle) = V(I) \cap V(f)$, i.e. that $f$ does not vanish on $V(I)$. Since $I$ is radical, the ideal of polynomials in $R$ vanishing on $V(I)$ is exactly $I$ (by the Nullstellensatz). By definition of zero-dimensional, the quotient $R/I$ is a finite-dimensional $\mathbb{C}$-vector space.</p>
<ul>
<li><p>Consider the $\mathbb{C}$-linear "multiplication by $f$" map $M_f: g \mapsto f \cdot g$ on $R/I$. Let $g$ be an eigenvector of $M_f$ with eigenvalue $\lambda$, so $f \cdot g = \lambda \cdot g$ or $(f-\lambda) g = 0$ holds in $R/I$, i.e. $(f-\lambda)g \in I$. The fact that $g \neq 0$ in $R/I$ means $g \not\in I$, hence there is a point $p \in V(I)$ with $g(p) \neq 0$. Now from $(f-\lambda)g \in I$ it follows that $\lambda = f(p)$ is the value of $f$ at the point $p$.</p></li>
<li><p>Conversely, if $P$ is the minimum polynomial of $M_f$ then $P(M_f) = 0$ implies $P(f) = 0$ in $R/I$ (by evaluating at $1 \in R/I$), hence $P(f) \in I$ and it follows that $0 = P(f)(p) = P(f(p))$ for all $p \in V(I)$, i.e. all values $f(p)$ for $p \in V(I)$ are zeros of $P$, and hence they are eigenvalues of $M_f$.</p></li>
</ul>
<p>So the values of $f$ on $V(I)$ are exactly the eigenvalues of $M_f$, and it suffices to show that $0$ is not one of them.</p>
<p>Given a Groebner basis $G$ of $I$, the multivariate polynomial division algorithm implies that the set of all monomials not divisible by any leading terms of $G$ is a basis of $R/I$ as a $\mathbb{C}$-vector space (called the normal basis). The matrix of $M_f$ with respect to this basis is easily computed:</p>
<pre><code>def multiplication_matrix(f, I):
G = I.groebner_basis()
f = f.reduce(G)
NB = I.normal_basis()
M = Matrix(f.base_ring(), len(NB))
for (i, m1) in enumerate(NB):
g = (f*m1).reduce(G)
M.set_column(i, [g.monomial_coefficient(m2) for m2 in NB])
return M
</code></pre>
<p>We apply it to the problem at hand:</p>
<pre><code>sage: f = 5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125
sage: I = R.ideal(C1 + C2)
sage: M_f = multiplication_matrix(f, I)
sage: M_f.det() != 0
True
</code></pre>
<p>Hence $f$ does not attain the value $0$ on $V(I)$, so $\varnothing = V(I + \langle f \rangle)$ and $\langle 1 \rangle = I + \langle f \rangle$.</p>
<p>Since $M_f$ is invertible, we can find a polynomial representative of $f^{-1} + I$:</p>
<pre><code>sage: NB = list(I.normal_basis())
sage: len(NB), NB.index(1)
(196, 195)
sage: finv_coeffs = M_f \ vector([0]*195 + [1])
sage: finv = sum(c*m for (c,m) in zip(finv_coeffs, NB))
sage: (f*finv).reduce(I)
1
</code></pre>
<p>Then we have $1 = f^{-1} \cdot f$ in $R/I$, i.e. $1 = f^{-1} \cdot f + g$ for some $g \in I$, and we can express $g$ in terms of generators by applying the division algorithm to $1 - f^{-1} \cdot f$ (or finding a lift of it to $I$):</p>
<pre><code>sage: lift_coeffs = [finv] + list((1 - finv*f).lift(I))
sage: sum(c*m for c,m in zip(lift_coeffs, [f] + I.gens()))
1
</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57248#post-id-57248Very good!
Then about the existence of solutions in char. $p>0$, in this case:
- first we need to consider independently the finite set $S$ of prime divisors of the (reduced) denominators in the Groebner basis of $I$,
- next, out of $S$, there is a solution in char. $p$ if $p$ divides the (reduced) numerator $n$ of $\det(M_f)$.
This "if" is not "iff" because both $41 \not \in S$ and does not divides $n$. In fact, the ideal $I$ is not zero-dimensional in char. $41$ but positive-dimensional (in fact 2-dim). I guess that $41$ is the only prime in this alternative case.Mon, 24 May 2021 04:59:30 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57248#post-id-57248Answer by mwageringel for <p>I have two (equivalent) sub-systems of polynomial equations of 10 variables and 12 equations, linked with one extra equations (so a total of 20 variable and 25 equations).</p>
<p>I can very quicky compute the Groebner basis of each sub-system. Then I can put them together, with the extra equation, and then I can quicky get the Groebner basis of the full system, which turns out to be trivial (total computation time: less than 3 minutes). See all the details in Appendix.</p>
<p>A direct computation of the Groebner basis of the full system is too long (I stopped after 1 day). The problem is that I need to compute the lift of 1 of the original full system (i.e. express 1 as a linear combination (with polynomial coefficients) of the original generators, by using the lift command <code>list(R.one().lift(Id))</code>), but it is also too long... </p>
<p>Now I expect that the application of a strategy involving the two sub-systems (as above) should exist to quickly get the lift of 1, but <strong>how?</strong> </p>
<p>I tried to compute this lift in two times, using the Groebner basis of the sub-systems, but it is also too long...</p>
<hr>
<h2>Appendix</h2>
<p>First sub-system (10 variables and 12 equations)</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
sage: A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
....: 5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
....: 25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
....: 5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
....: 5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
....: 5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
....: 5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
....: 25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
....: 5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
....: 25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
....: 5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
....: 5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
sage: Id=R.ideal(A1)
sage: %time G1=Id.groebner_basis()
CPU times: user 1.33 s, sys: 0 ns, total: 1.33 s
Wall time: 1.33 s
sage: C1=[g for g in G1]
</code></pre>
<p>Second (equivalent) sub-system: </p>
<pre><code>sage: R.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
sage: A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
....: 5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
....: 25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
....: 5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
....: 5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
....: 5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
....: 5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
....: 25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
....: 5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
....: 5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
....: 25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
....: 5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
sage: Id=R.ideal(A2)
sage: %time G2=Id.groebner_basis()
CPU times: user 2.08 s, sys: 0 ns, total: 2.08 s
Wall time: 2.08 s
sage: C2=[g for g in G2]
</code></pre>
<p>The full system:</p>
<pre><code>sage: R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
sage: C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
sage: Id=R.ideal(C)
sage: %time G=Id.groebner_basis()
CPU times: user 2min 29s, sys: 16 ms, total: 2min 29s
Wall time: 2min 29s
sage: G
[1]
</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?answer=57233#post-id-57233Note that computing the lift of 1 with respect to your `C` modulo a prime is fast – it takes only a few seconds. Doing this for many primes can allow to lift the result to ℚ. In this case, the result `F` (<s>see attachment</s> [uploaded here (31MB)](https://myshare.uni-osnabrueck.de/f/cfe7fa26c18d47c691c8/?dl=1), for now, as the file seems too large for an attachment) is huge, though, with coefficients as large as 4826 and 4818 digits for numerator and denominator.
<pre><code>sage: F, C = load("lift_one_result_F_C.sobj") # saved with Sage 9.3
sage: sum(f*c for f,c in zip(F, C)) # 1 minute
1</code></pre>
I am not sure if this solution would be unique. Maybe there is a smaller one.
To obtain the result with respect to `A1` and `A2`, you could try to lift `C1` and `C2` individually, but this might make the result even larger.
<hr>
Update
----
When computing the lift of 1 mod p, most of the polynomials in the result have the same monomials for their support, if p varies. This makes it very likely that the correct representative over ℚ has those same monomials as support and that the polynomials we computed mod p are the mod-p versions of the desired polynomials over ℚ.
This can fail for some primes which we need to ignore. In the function `likely_support`, we use a simple ad-hoc heuristic to find the exponents that probably occur in the ℚ-result: if a monomial occurs in 50% of the mod-p polynomials (i.e. modulo 50% of the primes or more), we expect that it should occur in the ℚ-result as well. The 50% threshold does not matter much. We could have picked 10% or 90% instead, I suppose.
We only keep the primes p for which all the polynomials in the lift of 1 mod p have the expected supports. Out of the 2800 primes we started with, 2358 remain.
Next, we use the Chinese remainder theorem `CRT_list` to lift the coefficients to ℤ modulo the product of the primes and then `rational_reconstruction` to obtain rational coefficients, both implemented in `lift_poly_crt`. For `rational_reconstruction` to work, the product of the primes must have about twice as many digits as the numerators and denominators in the result over ℚ (which I do not know how to determine in advance). In our case, we obtain a modulus with 10538 digits.
I have used `cached_function.precompute` in order to parallelize some computations because it is a bit less awkward to use than the `@parallel` decorator, but also because it allows to increase the number of primes until rational reconstruction works, without redoing all the previous slow computations.
The slowest part of the computation is lifting 1 modulo the 2800 primes, so this benefits a lot from parallelization. Each prime takes about 10 seconds, so this would take about 7.5 hours on a single thread.
<pre><code>R1.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
I1=R1.ideal(A1)
C1=list(I1.groebner_basis())
R2.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
I2=R2.ideal(A2)
C2=list(I2.groebner_basis())
R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
I=R.ideal(C)
########################################################################
@cached_function
def lift_one_modp(Ip):
Fp = Ip.ring()(1).lift(Ip.ideal())
print("completed p = %s" % Fp.ring().characteristic())
return Fp
@cached_function
def lift_poly_crt(fps):
ps = [fp.parent().characteristic() for fp in fps]
modulus = prod(ps)
Ds = [fp.dict() for fp in fps]
supp = set([e for Di in Ds for e in Di])
return R({e: CRT_list([Di[e].lift() if e in Di else 0 for Di in Ds],
list(ps)).rational_reconstruction(modulus) for e in supp})
def likely_support(fps):
supps = [set(fp.dict().keys()) for fp in fps]
supp_union = set.union(*supps)
supp = set()
half = len(supps) * 1/2 # exact amount does not matter
for e in supp_union:
k = 0
for s in supps:
if e in s:
k += 1
if k > half: # if exponent vector is in more than half the supports, we will probably want to keep it
supp.add(e)
break
return supp
ps = sorted(RecursivelyEnumeratedSet([next_prime(16000)], lambda p: [next_prime(p)])[:2800]) # many primes
Rps = [R.change_ring(GF(p)) for p in ps]
Ips = [I.change_ring(Rp).gens() for Rp in Rps]
CPUs = 16 # choose based on hardware
lift_one_modp.precompute(Ips, num_processes=CPUs)
assert all(type(r) != str for r in lift_one_modp.cache.values()) # otherwise, there was an exception
Fps = [lift_one_modp(Ip) for Ip in Ips]
supps = [likely_support(fps) for fps in zip(*Fps)]
# keep indices of primes for which all polynomials have the expected support
idxs = [j for j, Fp in enumerate(Fps) if all(set(fp.dict()) == supp for fp, supp in zip(Fp, supps))]
# restrict the results to those primes
ps1 = tuple(ps[j] for j in idxs)
Fps1 = [Fps[j] for j in idxs]
fpss = [tuple(fps) for fps in zip(*Fps1)]
modulus1 = prod(ps1)
floor(log(modulus1, 10)) # 10538 (must be large enough for rational reconstruction: this needs about twice as many digits as the numerators or denominators in the result have)
lift_poly_crt.precompute([(fps,) for fps in fpss], num_processes=CPUs)
assert all(type(r) != str for r in lift_poly_crt.cache.values()) # otherwise, there was an exception (if so, remove it from cache by lift_poly_crt.cache.clear())
F = [lift_poly_crt(fps) for fps in fpss] # the final result</code></pre>Sat, 22 May 2021 23:52:48 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?answer=57233#post-id-57233Comment by Sébastien Palcoux for <p>Note that computing the lift of 1 with respect to your <code>C</code> modulo a prime is fast – it takes only a few seconds. Doing this for many primes can allow to lift the result to ℚ. In this case, the result <code>F</code> (<s>see attachment</s> <a href="https://myshare.uni-osnabrueck.de/f/cfe7fa26c18d47c691c8/?dl=1">uploaded here (31MB)</a>, for now, as the file seems too large for an attachment) is huge, though, with coefficients as large as 4826 and 4818 digits for numerator and denominator.</p>
<pre><code>sage: F, C = load("lift_one_result_F_C.sobj") # saved with Sage 9.3
sage: sum(f*c for f,c in zip(F, C)) # 1 minute
1</code></pre>
<p>I am not sure if this solution would be unique. Maybe there is a smaller one.</p>
<p>To obtain the result with respect to <code>A1</code> and <code>A2</code>, you could try to lift <code>C1</code> and <code>C2</code> individually, but this might make the result even larger.</p>
<hr>
<h2>Update</h2>
<p>When computing the lift of 1 mod p, most of the polynomials in the result have the same monomials for their support, if p varies. This makes it very likely that the correct representative over ℚ has those same monomials as support and that the polynomials we computed mod p are the mod-p versions of the desired polynomials over ℚ.</p>
<p>This can fail for some primes which we need to ignore. In the function <code>likely_support</code>, we use a simple ad-hoc heuristic to find the exponents that probably occur in the ℚ-result: if a monomial occurs in 50% of the mod-p polynomials (i.e. modulo 50% of the primes or more), we expect that it should occur in the ℚ-result as well. The 50% threshold does not matter much. We could have picked 10% or 90% instead, I suppose.</p>
<p>We only keep the primes p for which all the polynomials in the lift of 1 mod p have the expected supports. Out of the 2800 primes we started with, 2358 remain.</p>
<p>Next, we use the Chinese remainder theorem <code>CRT_list</code> to lift the coefficients to ℤ modulo the product of the primes and then <code>rational_reconstruction</code> to obtain rational coefficients, both implemented in <code>lift_poly_crt</code>. For <code>rational_reconstruction</code> to work, the product of the primes must have about twice as many digits as the numerators and denominators in the result over ℚ (which I do not know how to determine in advance). In our case, we obtain a modulus with 10538 digits.</p>
<p>I have used <code>cached_function.precompute</code> in order to parallelize some computations because it is a bit less awkward to use than the <code>@parallel</code> decorator, but also because it allows to increase the number of primes until rational reconstruction works, without redoing all the previous slow computations.</p>
<p>The slowest part of the computation is lifting 1 modulo the 2800 primes, so this benefits a lot from parallelization. Each prime takes about 10 seconds, so this would take about 7.5 hours on a single thread.</p>
<pre><code>R1.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
I1=R1.ideal(A1)
C1=list(I1.groebner_basis())
R2.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
I2=R2.ideal(A2)
C2=list(I2.groebner_basis())
R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
I=R.ideal(C)
########################################################################
@cached_function
def lift_one_modp(Ip):
Fp = Ip.ring()(1).lift(Ip.ideal())
print("completed p = %s" % Fp.ring().characteristic())
return Fp
@cached_function
def lift_poly_crt(fps):
ps = [fp.parent().characteristic() for fp in fps]
modulus = prod(ps)
Ds = [fp.dict() for fp in fps]
supp = set([e for Di in Ds for e in Di])
return R({e: CRT_list([Di[e].lift() if e in Di else 0 for Di in Ds],
list(ps)).rational_reconstruction(modulus) for e in supp})
def likely_support(fps):
supps = [set(fp.dict().keys()) for fp in fps]
supp_union = set.union(*supps)
supp = set()
half = len(supps) * 1/2 # exact amount does not matter
for e in supp_union:
k = 0
for s in supps:
if e in s:
k += 1
if k > half: # if exponent vector is in more than half the supports, we will probably want to keep it
supp.add(e)
break
return supp
ps = sorted(RecursivelyEnumeratedSet([next_prime(16000)], lambda p: [next_prime(p)])[:2800]) # many primes
Rps = [R.change_ring(GF(p)) for p in ps]
Ips = [I.change_ring(Rp).gens() for Rp in Rps]
CPUs = 16 # choose based on hardware
lift_one_modp.precompute(Ips, num_processes=CPUs)
assert all(type(r) != str for r in lift_one_modp.cache.values()) # otherwise, there was an exception
Fps = [lift_one_modp(Ip) for Ip in Ips]
supps = [likely_support(fps) for fps in zip(*Fps)]
# keep indices of primes for which all polynomials have the expected support
idxs = [j for j, Fp in enumerate(Fps) if all(set(fp.dict()) == supp for fp, supp in zip(Fp, supps))]
# restrict the results to those primes
ps1 = tuple(ps[j] for j in idxs)
Fps1 = [Fps[j] for j in idxs]
fpss = [tuple(fps) for fps in zip(*Fps1)]
modulus1 = prod(ps1)
floor(log(modulus1, 10)) # 10538 (must be large enough for rational reconstruction: this needs about twice as many digits as the numerators or denominators in the result have)
lift_poly_crt.precompute([(fps,) for fps in fpss], num_processes=CPUs)
assert all(type(r) != str for r in lift_poly_crt.cache.values()) # otherwise, there was an exception (if so, remove it from cache by lift_poly_crt.cache.clear())
F = [lift_poly_crt(fps) for fps in fpss] # the final result</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57935#post-id-57935@mwageringel: see [this new post](https://ask.sagemath.org/question/57934/polynomial-system-without-solution-in-char-0-classification-of-char-p-with-solution/) dedicated to the problem in my previous comment.Fri, 09 Jul 2021 11:29:33 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57935#post-id-57935Comment by Sébastien Palcoux for <p>Note that computing the lift of 1 with respect to your <code>C</code> modulo a prime is fast – it takes only a few seconds. Doing this for many primes can allow to lift the result to ℚ. In this case, the result <code>F</code> (<s>see attachment</s> <a href="https://myshare.uni-osnabrueck.de/f/cfe7fa26c18d47c691c8/?dl=1">uploaded here (31MB)</a>, for now, as the file seems too large for an attachment) is huge, though, with coefficients as large as 4826 and 4818 digits for numerator and denominator.</p>
<pre><code>sage: F, C = load("lift_one_result_F_C.sobj") # saved with Sage 9.3
sage: sum(f*c for f,c in zip(F, C)) # 1 minute
1</code></pre>
<p>I am not sure if this solution would be unique. Maybe there is a smaller one.</p>
<p>To obtain the result with respect to <code>A1</code> and <code>A2</code>, you could try to lift <code>C1</code> and <code>C2</code> individually, but this might make the result even larger.</p>
<hr>
<h2>Update</h2>
<p>When computing the lift of 1 mod p, most of the polynomials in the result have the same monomials for their support, if p varies. This makes it very likely that the correct representative over ℚ has those same monomials as support and that the polynomials we computed mod p are the mod-p versions of the desired polynomials over ℚ.</p>
<p>This can fail for some primes which we need to ignore. In the function <code>likely_support</code>, we use a simple ad-hoc heuristic to find the exponents that probably occur in the ℚ-result: if a monomial occurs in 50% of the mod-p polynomials (i.e. modulo 50% of the primes or more), we expect that it should occur in the ℚ-result as well. The 50% threshold does not matter much. We could have picked 10% or 90% instead, I suppose.</p>
<p>We only keep the primes p for which all the polynomials in the lift of 1 mod p have the expected supports. Out of the 2800 primes we started with, 2358 remain.</p>
<p>Next, we use the Chinese remainder theorem <code>CRT_list</code> to lift the coefficients to ℤ modulo the product of the primes and then <code>rational_reconstruction</code> to obtain rational coefficients, both implemented in <code>lift_poly_crt</code>. For <code>rational_reconstruction</code> to work, the product of the primes must have about twice as many digits as the numerators and denominators in the result over ℚ (which I do not know how to determine in advance). In our case, we obtain a modulus with 10538 digits.</p>
<p>I have used <code>cached_function.precompute</code> in order to parallelize some computations because it is a bit less awkward to use than the <code>@parallel</code> decorator, but also because it allows to increase the number of primes until rational reconstruction works, without redoing all the previous slow computations.</p>
<p>The slowest part of the computation is lifting 1 modulo the 2800 primes, so this benefits a lot from parallelization. Each prime takes about 10 seconds, so this would take about 7.5 hours on a single thread.</p>
<pre><code>R1.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
I1=R1.ideal(A1)
C1=list(I1.groebner_basis())
R2.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
I2=R2.ideal(A2)
C2=list(I2.groebner_basis())
R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
I=R.ideal(C)
########################################################################
@cached_function
def lift_one_modp(Ip):
Fp = Ip.ring()(1).lift(Ip.ideal())
print("completed p = %s" % Fp.ring().characteristic())
return Fp
@cached_function
def lift_poly_crt(fps):
ps = [fp.parent().characteristic() for fp in fps]
modulus = prod(ps)
Ds = [fp.dict() for fp in fps]
supp = set([e for Di in Ds for e in Di])
return R({e: CRT_list([Di[e].lift() if e in Di else 0 for Di in Ds],
list(ps)).rational_reconstruction(modulus) for e in supp})
def likely_support(fps):
supps = [set(fp.dict().keys()) for fp in fps]
supp_union = set.union(*supps)
supp = set()
half = len(supps) * 1/2 # exact amount does not matter
for e in supp_union:
k = 0
for s in supps:
if e in s:
k += 1
if k > half: # if exponent vector is in more than half the supports, we will probably want to keep it
supp.add(e)
break
return supp
ps = sorted(RecursivelyEnumeratedSet([next_prime(16000)], lambda p: [next_prime(p)])[:2800]) # many primes
Rps = [R.change_ring(GF(p)) for p in ps]
Ips = [I.change_ring(Rp).gens() for Rp in Rps]
CPUs = 16 # choose based on hardware
lift_one_modp.precompute(Ips, num_processes=CPUs)
assert all(type(r) != str for r in lift_one_modp.cache.values()) # otherwise, there was an exception
Fps = [lift_one_modp(Ip) for Ip in Ips]
supps = [likely_support(fps) for fps in zip(*Fps)]
# keep indices of primes for which all polynomials have the expected support
idxs = [j for j, Fp in enumerate(Fps) if all(set(fp.dict()) == supp for fp, supp in zip(Fp, supps))]
# restrict the results to those primes
ps1 = tuple(ps[j] for j in idxs)
Fps1 = [Fps[j] for j in idxs]
fpss = [tuple(fps) for fps in zip(*Fps1)]
modulus1 = prod(ps1)
floor(log(modulus1, 10)) # 10538 (must be large enough for rational reconstruction: this needs about twice as many digits as the numerators or denominators in the result have)
lift_poly_crt.precompute([(fps,) for fps in fpss], num_processes=CPUs)
assert all(type(r) != str for r in lift_poly_crt.cache.values()) # otherwise, there was an exception (if so, remove it from cache by lift_poly_crt.cache.clear())
F = [lift_poly_crt(fps) for fps in fpss] # the final result</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57887#post-id-57887@mwageringel: in fact, what I really need is the full list of prime numbers p such that the system A1+A2+extra admits a non-trivial Groebner basis over GF(p). Some of them are the prime factors of the numerator of the determinant given by rburing, but not all (for example p=41). To get all of them, my strategy was to lift C1 with respect to A1, and to test all the prime divisors of the denominators there.Tue, 06 Jul 2021 15:31:04 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57887#post-id-57887Comment by mwageringel for <p>Note that computing the lift of 1 with respect to your <code>C</code> modulo a prime is fast – it takes only a few seconds. Doing this for many primes can allow to lift the result to ℚ. In this case, the result <code>F</code> (<s>see attachment</s> <a href="https://myshare.uni-osnabrueck.de/f/cfe7fa26c18d47c691c8/?dl=1">uploaded here (31MB)</a>, for now, as the file seems too large for an attachment) is huge, though, with coefficients as large as 4826 and 4818 digits for numerator and denominator.</p>
<pre><code>sage: F, C = load("lift_one_result_F_C.sobj") # saved with Sage 9.3
sage: sum(f*c for f,c in zip(F, C)) # 1 minute
1</code></pre>
<p>I am not sure if this solution would be unique. Maybe there is a smaller one.</p>
<p>To obtain the result with respect to <code>A1</code> and <code>A2</code>, you could try to lift <code>C1</code> and <code>C2</code> individually, but this might make the result even larger.</p>
<hr>
<h2>Update</h2>
<p>When computing the lift of 1 mod p, most of the polynomials in the result have the same monomials for their support, if p varies. This makes it very likely that the correct representative over ℚ has those same monomials as support and that the polynomials we computed mod p are the mod-p versions of the desired polynomials over ℚ.</p>
<p>This can fail for some primes which we need to ignore. In the function <code>likely_support</code>, we use a simple ad-hoc heuristic to find the exponents that probably occur in the ℚ-result: if a monomial occurs in 50% of the mod-p polynomials (i.e. modulo 50% of the primes or more), we expect that it should occur in the ℚ-result as well. The 50% threshold does not matter much. We could have picked 10% or 90% instead, I suppose.</p>
<p>We only keep the primes p for which all the polynomials in the lift of 1 mod p have the expected supports. Out of the 2800 primes we started with, 2358 remain.</p>
<p>Next, we use the Chinese remainder theorem <code>CRT_list</code> to lift the coefficients to ℤ modulo the product of the primes and then <code>rational_reconstruction</code> to obtain rational coefficients, both implemented in <code>lift_poly_crt</code>. For <code>rational_reconstruction</code> to work, the product of the primes must have about twice as many digits as the numerators and denominators in the result over ℚ (which I do not know how to determine in advance). In our case, we obtain a modulus with 10538 digits.</p>
<p>I have used <code>cached_function.precompute</code> in order to parallelize some computations because it is a bit less awkward to use than the <code>@parallel</code> decorator, but also because it allows to increase the number of primes until rational reconstruction works, without redoing all the previous slow computations.</p>
<p>The slowest part of the computation is lifting 1 modulo the 2800 primes, so this benefits a lot from parallelization. Each prime takes about 10 seconds, so this would take about 7.5 hours on a single thread.</p>
<pre><code>R1.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
I1=R1.ideal(A1)
C1=list(I1.groebner_basis())
R2.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
I2=R2.ideal(A2)
C2=list(I2.groebner_basis())
R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
I=R.ideal(C)
########################################################################
@cached_function
def lift_one_modp(Ip):
Fp = Ip.ring()(1).lift(Ip.ideal())
print("completed p = %s" % Fp.ring().characteristic())
return Fp
@cached_function
def lift_poly_crt(fps):
ps = [fp.parent().characteristic() for fp in fps]
modulus = prod(ps)
Ds = [fp.dict() for fp in fps]
supp = set([e for Di in Ds for e in Di])
return R({e: CRT_list([Di[e].lift() if e in Di else 0 for Di in Ds],
list(ps)).rational_reconstruction(modulus) for e in supp})
def likely_support(fps):
supps = [set(fp.dict().keys()) for fp in fps]
supp_union = set.union(*supps)
supp = set()
half = len(supps) * 1/2 # exact amount does not matter
for e in supp_union:
k = 0
for s in supps:
if e in s:
k += 1
if k > half: # if exponent vector is in more than half the supports, we will probably want to keep it
supp.add(e)
break
return supp
ps = sorted(RecursivelyEnumeratedSet([next_prime(16000)], lambda p: [next_prime(p)])[:2800]) # many primes
Rps = [R.change_ring(GF(p)) for p in ps]
Ips = [I.change_ring(Rp).gens() for Rp in Rps]
CPUs = 16 # choose based on hardware
lift_one_modp.precompute(Ips, num_processes=CPUs)
assert all(type(r) != str for r in lift_one_modp.cache.values()) # otherwise, there was an exception
Fps = [lift_one_modp(Ip) for Ip in Ips]
supps = [likely_support(fps) for fps in zip(*Fps)]
# keep indices of primes for which all polynomials have the expected support
idxs = [j for j, Fp in enumerate(Fps) if all(set(fp.dict()) == supp for fp, supp in zip(Fp, supps))]
# restrict the results to those primes
ps1 = tuple(ps[j] for j in idxs)
Fps1 = [Fps[j] for j in idxs]
fpss = [tuple(fps) for fps in zip(*Fps1)]
modulus1 = prod(ps1)
floor(log(modulus1, 10)) # 10538 (must be large enough for rational reconstruction: this needs about twice as many digits as the numerators or denominators in the result have)
lift_poly_crt.precompute([(fps,) for fps in fpss], num_processes=CPUs)
assert all(type(r) != str for r in lift_poly_crt.cache.values()) # otherwise, there was an exception (if so, remove it from cache by lift_poly_crt.cache.clear())
F = [lift_poly_crt(fps) for fps in fpss] # the final result</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57854#post-id-57854@Sébastien Palcoux Instead of writing `1` as a linear combination of the generators of `C`, you would want to write each generator of `C1` as a linear combination of the generators of `A1`. So `A1`, `I1`, `R1` take the roles of `C`, `I`, `R` respectively, and in `lift_one_modp` replace the `1` by, say, `C1[0]`, the first generator of `C1`. This allows you to express `C1[0]` as a linear combination of `A1`, which then needs to be repeated for all the other generators of `C1`.
If I recall correctly, I have tried this computation for a single prime, which was already very slow and the result was very large, so repeating that for lots of primes in order to lift the result to ℚ may not be viable.Fri, 02 Jul 2021 20:42:44 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57854#post-id-57854Comment by Sébastien Palcoux for <p>Note that computing the lift of 1 with respect to your <code>C</code> modulo a prime is fast – it takes only a few seconds. Doing this for many primes can allow to lift the result to ℚ. In this case, the result <code>F</code> (<s>see attachment</s> <a href="https://myshare.uni-osnabrueck.de/f/cfe7fa26c18d47c691c8/?dl=1">uploaded here (31MB)</a>, for now, as the file seems too large for an attachment) is huge, though, with coefficients as large as 4826 and 4818 digits for numerator and denominator.</p>
<pre><code>sage: F, C = load("lift_one_result_F_C.sobj") # saved with Sage 9.3
sage: sum(f*c for f,c in zip(F, C)) # 1 minute
1</code></pre>
<p>I am not sure if this solution would be unique. Maybe there is a smaller one.</p>
<p>To obtain the result with respect to <code>A1</code> and <code>A2</code>, you could try to lift <code>C1</code> and <code>C2</code> individually, but this might make the result even larger.</p>
<hr>
<h2>Update</h2>
<p>When computing the lift of 1 mod p, most of the polynomials in the result have the same monomials for their support, if p varies. This makes it very likely that the correct representative over ℚ has those same monomials as support and that the polynomials we computed mod p are the mod-p versions of the desired polynomials over ℚ.</p>
<p>This can fail for some primes which we need to ignore. In the function <code>likely_support</code>, we use a simple ad-hoc heuristic to find the exponents that probably occur in the ℚ-result: if a monomial occurs in 50% of the mod-p polynomials (i.e. modulo 50% of the primes or more), we expect that it should occur in the ℚ-result as well. The 50% threshold does not matter much. We could have picked 10% or 90% instead, I suppose.</p>
<p>We only keep the primes p for which all the polynomials in the lift of 1 mod p have the expected supports. Out of the 2800 primes we started with, 2358 remain.</p>
<p>Next, we use the Chinese remainder theorem <code>CRT_list</code> to lift the coefficients to ℤ modulo the product of the primes and then <code>rational_reconstruction</code> to obtain rational coefficients, both implemented in <code>lift_poly_crt</code>. For <code>rational_reconstruction</code> to work, the product of the primes must have about twice as many digits as the numerators and denominators in the result over ℚ (which I do not know how to determine in advance). In our case, we obtain a modulus with 10538 digits.</p>
<p>I have used <code>cached_function.precompute</code> in order to parallelize some computations because it is a bit less awkward to use than the <code>@parallel</code> decorator, but also because it allows to increase the number of primes until rational reconstruction works, without redoing all the previous slow computations.</p>
<p>The slowest part of the computation is lifting 1 modulo the 2800 primes, so this benefits a lot from parallelization. Each prime takes about 10 seconds, so this would take about 7.5 hours on a single thread.</p>
<pre><code>R1.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
I1=R1.ideal(A1)
C1=list(I1.groebner_basis())
R2.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
I2=R2.ideal(A2)
C2=list(I2.groebner_basis())
R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
I=R.ideal(C)
########################################################################
@cached_function
def lift_one_modp(Ip):
Fp = Ip.ring()(1).lift(Ip.ideal())
print("completed p = %s" % Fp.ring().characteristic())
return Fp
@cached_function
def lift_poly_crt(fps):
ps = [fp.parent().characteristic() for fp in fps]
modulus = prod(ps)
Ds = [fp.dict() for fp in fps]
supp = set([e for Di in Ds for e in Di])
return R({e: CRT_list([Di[e].lift() if e in Di else 0 for Di in Ds],
list(ps)).rational_reconstruction(modulus) for e in supp})
def likely_support(fps):
supps = [set(fp.dict().keys()) for fp in fps]
supp_union = set.union(*supps)
supp = set()
half = len(supps) * 1/2 # exact amount does not matter
for e in supp_union:
k = 0
for s in supps:
if e in s:
k += 1
if k > half: # if exponent vector is in more than half the supports, we will probably want to keep it
supp.add(e)
break
return supp
ps = sorted(RecursivelyEnumeratedSet([next_prime(16000)], lambda p: [next_prime(p)])[:2800]) # many primes
Rps = [R.change_ring(GF(p)) for p in ps]
Ips = [I.change_ring(Rp).gens() for Rp in Rps]
CPUs = 16 # choose based on hardware
lift_one_modp.precompute(Ips, num_processes=CPUs)
assert all(type(r) != str for r in lift_one_modp.cache.values()) # otherwise, there was an exception
Fps = [lift_one_modp(Ip) for Ip in Ips]
supps = [likely_support(fps) for fps in zip(*Fps)]
# keep indices of primes for which all polynomials have the expected support
idxs = [j for j, Fp in enumerate(Fps) if all(set(fp.dict()) == supp for fp, supp in zip(Fp, supps))]
# restrict the results to those primes
ps1 = tuple(ps[j] for j in idxs)
Fps1 = [Fps[j] for j in idxs]
fpss = [tuple(fps) for fps in zip(*Fps1)]
modulus1 = prod(ps1)
floor(log(modulus1, 10)) # 10538 (must be large enough for rational reconstruction: this needs about twice as many digits as the numerators or denominators in the result have)
lift_poly_crt.precompute([(fps,) for fps in fpss], num_processes=CPUs)
assert all(type(r) != str for r in lift_poly_crt.cache.values()) # otherwise, there was an exception (if so, remove it from cache by lift_poly_crt.cache.clear())
F = [lift_poly_crt(fps) for fps in fpss] # the final result</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57849#post-id-57849@mwageringel How can I use your code to lift $C_1$ with respect to $A_1$?Fri, 02 Jul 2021 11:24:53 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57849#post-id-57849Comment by mwageringel for <p>Note that computing the lift of 1 with respect to your <code>C</code> modulo a prime is fast – it takes only a few seconds. Doing this for many primes can allow to lift the result to ℚ. In this case, the result <code>F</code> (<s>see attachment</s> <a href="https://myshare.uni-osnabrueck.de/f/cfe7fa26c18d47c691c8/?dl=1">uploaded here (31MB)</a>, for now, as the file seems too large for an attachment) is huge, though, with coefficients as large as 4826 and 4818 digits for numerator and denominator.</p>
<pre><code>sage: F, C = load("lift_one_result_F_C.sobj") # saved with Sage 9.3
sage: sum(f*c for f,c in zip(F, C)) # 1 minute
1</code></pre>
<p>I am not sure if this solution would be unique. Maybe there is a smaller one.</p>
<p>To obtain the result with respect to <code>A1</code> and <code>A2</code>, you could try to lift <code>C1</code> and <code>C2</code> individually, but this might make the result even larger.</p>
<hr>
<h2>Update</h2>
<p>When computing the lift of 1 mod p, most of the polynomials in the result have the same monomials for their support, if p varies. This makes it very likely that the correct representative over ℚ has those same monomials as support and that the polynomials we computed mod p are the mod-p versions of the desired polynomials over ℚ.</p>
<p>This can fail for some primes which we need to ignore. In the function <code>likely_support</code>, we use a simple ad-hoc heuristic to find the exponents that probably occur in the ℚ-result: if a monomial occurs in 50% of the mod-p polynomials (i.e. modulo 50% of the primes or more), we expect that it should occur in the ℚ-result as well. The 50% threshold does not matter much. We could have picked 10% or 90% instead, I suppose.</p>
<p>We only keep the primes p for which all the polynomials in the lift of 1 mod p have the expected supports. Out of the 2800 primes we started with, 2358 remain.</p>
<p>Next, we use the Chinese remainder theorem <code>CRT_list</code> to lift the coefficients to ℤ modulo the product of the primes and then <code>rational_reconstruction</code> to obtain rational coefficients, both implemented in <code>lift_poly_crt</code>. For <code>rational_reconstruction</code> to work, the product of the primes must have about twice as many digits as the numerators and denominators in the result over ℚ (which I do not know how to determine in advance). In our case, we obtain a modulus with 10538 digits.</p>
<p>I have used <code>cached_function.precompute</code> in order to parallelize some computations because it is a bit less awkward to use than the <code>@parallel</code> decorator, but also because it allows to increase the number of primes until rational reconstruction works, without redoing all the previous slow computations.</p>
<p>The slowest part of the computation is lifting 1 modulo the 2800 primes, so this benefits a lot from parallelization. Each prime takes about 10 seconds, so this would take about 7.5 hours on a single thread.</p>
<pre><code>R1.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
I1=R1.ideal(A1)
C1=list(I1.groebner_basis())
R2.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
I2=R2.ideal(A2)
C2=list(I2.groebner_basis())
R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
I=R.ideal(C)
########################################################################
@cached_function
def lift_one_modp(Ip):
Fp = Ip.ring()(1).lift(Ip.ideal())
print("completed p = %s" % Fp.ring().characteristic())
return Fp
@cached_function
def lift_poly_crt(fps):
ps = [fp.parent().characteristic() for fp in fps]
modulus = prod(ps)
Ds = [fp.dict() for fp in fps]
supp = set([e for Di in Ds for e in Di])
return R({e: CRT_list([Di[e].lift() if e in Di else 0 for Di in Ds],
list(ps)).rational_reconstruction(modulus) for e in supp})
def likely_support(fps):
supps = [set(fp.dict().keys()) for fp in fps]
supp_union = set.union(*supps)
supp = set()
half = len(supps) * 1/2 # exact amount does not matter
for e in supp_union:
k = 0
for s in supps:
if e in s:
k += 1
if k > half: # if exponent vector is in more than half the supports, we will probably want to keep it
supp.add(e)
break
return supp
ps = sorted(RecursivelyEnumeratedSet([next_prime(16000)], lambda p: [next_prime(p)])[:2800]) # many primes
Rps = [R.change_ring(GF(p)) for p in ps]
Ips = [I.change_ring(Rp).gens() for Rp in Rps]
CPUs = 16 # choose based on hardware
lift_one_modp.precompute(Ips, num_processes=CPUs)
assert all(type(r) != str for r in lift_one_modp.cache.values()) # otherwise, there was an exception
Fps = [lift_one_modp(Ip) for Ip in Ips]
supps = [likely_support(fps) for fps in zip(*Fps)]
# keep indices of primes for which all polynomials have the expected support
idxs = [j for j, Fp in enumerate(Fps) if all(set(fp.dict()) == supp for fp, supp in zip(Fp, supps))]
# restrict the results to those primes
ps1 = tuple(ps[j] for j in idxs)
Fps1 = [Fps[j] for j in idxs]
fpss = [tuple(fps) for fps in zip(*Fps1)]
modulus1 = prod(ps1)
floor(log(modulus1, 10)) # 10538 (must be large enough for rational reconstruction: this needs about twice as many digits as the numerators or denominators in the result have)
lift_poly_crt.precompute([(fps,) for fps in fpss], num_processes=CPUs)
assert all(type(r) != str for r in lift_poly_crt.cache.values()) # otherwise, there was an exception (if so, remove it from cache by lift_poly_crt.cache.clear())
F = [lift_poly_crt(fps) for fps in fpss] # the final result</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57253#post-id-57253The trouble is that lifting $1$ with respect to $A_1$ and $A_2$ is too slow in general, even modulo $p$ it does not finish within a couple of minutes.
I have also tried to lift the elements of $C_1,C_2$ with respect to $A_1,A_2$ modulo $p$ – this finishes quickly, so doing this for many $p$ should allow lifting to ℚ again. However, if you compose this result with the $F$ from my answer, you get very large polynomials with about 500k terms. (Using `lift_coeffs` from @rburing's answer instead of $F$ would be a bit better, of course.)Mon, 24 May 2021 11:19:02 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57253#post-id-57253Comment by Sébastien Palcoux for <p>Note that computing the lift of 1 with respect to your <code>C</code> modulo a prime is fast – it takes only a few seconds. Doing this for many primes can allow to lift the result to ℚ. In this case, the result <code>F</code> (<s>see attachment</s> <a href="https://myshare.uni-osnabrueck.de/f/cfe7fa26c18d47c691c8/?dl=1">uploaded here (31MB)</a>, for now, as the file seems too large for an attachment) is huge, though, with coefficients as large as 4826 and 4818 digits for numerator and denominator.</p>
<pre><code>sage: F, C = load("lift_one_result_F_C.sobj") # saved with Sage 9.3
sage: sum(f*c for f,c in zip(F, C)) # 1 minute
1</code></pre>
<p>I am not sure if this solution would be unique. Maybe there is a smaller one.</p>
<p>To obtain the result with respect to <code>A1</code> and <code>A2</code>, you could try to lift <code>C1</code> and <code>C2</code> individually, but this might make the result even larger.</p>
<hr>
<h2>Update</h2>
<p>When computing the lift of 1 mod p, most of the polynomials in the result have the same monomials for their support, if p varies. This makes it very likely that the correct representative over ℚ has those same monomials as support and that the polynomials we computed mod p are the mod-p versions of the desired polynomials over ℚ.</p>
<p>This can fail for some primes which we need to ignore. In the function <code>likely_support</code>, we use a simple ad-hoc heuristic to find the exponents that probably occur in the ℚ-result: if a monomial occurs in 50% of the mod-p polynomials (i.e. modulo 50% of the primes or more), we expect that it should occur in the ℚ-result as well. The 50% threshold does not matter much. We could have picked 10% or 90% instead, I suppose.</p>
<p>We only keep the primes p for which all the polynomials in the lift of 1 mod p have the expected supports. Out of the 2800 primes we started with, 2358 remain.</p>
<p>Next, we use the Chinese remainder theorem <code>CRT_list</code> to lift the coefficients to ℤ modulo the product of the primes and then <code>rational_reconstruction</code> to obtain rational coefficients, both implemented in <code>lift_poly_crt</code>. For <code>rational_reconstruction</code> to work, the product of the primes must have about twice as many digits as the numerators and denominators in the result over ℚ (which I do not know how to determine in advance). In our case, we obtain a modulus with 10538 digits.</p>
<p>I have used <code>cached_function.precompute</code> in order to parallelize some computations because it is a bit less awkward to use than the <code>@parallel</code> decorator, but also because it allows to increase the number of primes until rational reconstruction works, without redoing all the previous slow computations.</p>
<p>The slowest part of the computation is lifting 1 modulo the 2800 primes, so this benefits a lot from parallelization. Each prime takes about 10 seconds, so this would take about 7.5 hours on a single thread.</p>
<pre><code>R1.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
I1=R1.ideal(A1)
C1=list(I1.groebner_basis())
R2.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
I2=R2.ideal(A2)
C2=list(I2.groebner_basis())
R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
I=R.ideal(C)
########################################################################
@cached_function
def lift_one_modp(Ip):
Fp = Ip.ring()(1).lift(Ip.ideal())
print("completed p = %s" % Fp.ring().characteristic())
return Fp
@cached_function
def lift_poly_crt(fps):
ps = [fp.parent().characteristic() for fp in fps]
modulus = prod(ps)
Ds = [fp.dict() for fp in fps]
supp = set([e for Di in Ds for e in Di])
return R({e: CRT_list([Di[e].lift() if e in Di else 0 for Di in Ds],
list(ps)).rational_reconstruction(modulus) for e in supp})
def likely_support(fps):
supps = [set(fp.dict().keys()) for fp in fps]
supp_union = set.union(*supps)
supp = set()
half = len(supps) * 1/2 # exact amount does not matter
for e in supp_union:
k = 0
for s in supps:
if e in s:
k += 1
if k > half: # if exponent vector is in more than half the supports, we will probably want to keep it
supp.add(e)
break
return supp
ps = sorted(RecursivelyEnumeratedSet([next_prime(16000)], lambda p: [next_prime(p)])[:2800]) # many primes
Rps = [R.change_ring(GF(p)) for p in ps]
Ips = [I.change_ring(Rp).gens() for Rp in Rps]
CPUs = 16 # choose based on hardware
lift_one_modp.precompute(Ips, num_processes=CPUs)
assert all(type(r) != str for r in lift_one_modp.cache.values()) # otherwise, there was an exception
Fps = [lift_one_modp(Ip) for Ip in Ips]
supps = [likely_support(fps) for fps in zip(*Fps)]
# keep indices of primes for which all polynomials have the expected support
idxs = [j for j, Fp in enumerate(Fps) if all(set(fp.dict()) == supp for fp, supp in zip(Fp, supps))]
# restrict the results to those primes
ps1 = tuple(ps[j] for j in idxs)
Fps1 = [Fps[j] for j in idxs]
fpss = [tuple(fps) for fps in zip(*Fps1)]
modulus1 = prod(ps1)
floor(log(modulus1, 10)) # 10538 (must be large enough for rational reconstruction: this needs about twice as many digits as the numerators or denominators in the result have)
lift_poly_crt.precompute([(fps,) for fps in fpss], num_processes=CPUs)
assert all(type(r) != str for r in lift_poly_crt.cache.values()) # otherwise, there was an exception (if so, remove it from cache by lift_poly_crt.cache.clear())
F = [lift_poly_crt(fps) for fps in fpss] # the final result</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57251#post-id-57251With a bit of luck, the lift of $1$ with respect to $A_1+A_2+$ extra could be easier...Mon, 24 May 2021 10:54:42 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57251#post-id-57251Comment by Sébastien Palcoux for <p>Note that computing the lift of 1 with respect to your <code>C</code> modulo a prime is fast – it takes only a few seconds. Doing this for many primes can allow to lift the result to ℚ. In this case, the result <code>F</code> (<s>see attachment</s> <a href="https://myshare.uni-osnabrueck.de/f/cfe7fa26c18d47c691c8/?dl=1">uploaded here (31MB)</a>, for now, as the file seems too large for an attachment) is huge, though, with coefficients as large as 4826 and 4818 digits for numerator and denominator.</p>
<pre><code>sage: F, C = load("lift_one_result_F_C.sobj") # saved with Sage 9.3
sage: sum(f*c for f,c in zip(F, C)) # 1 minute
1</code></pre>
<p>I am not sure if this solution would be unique. Maybe there is a smaller one.</p>
<p>To obtain the result with respect to <code>A1</code> and <code>A2</code>, you could try to lift <code>C1</code> and <code>C2</code> individually, but this might make the result even larger.</p>
<hr>
<h2>Update</h2>
<p>When computing the lift of 1 mod p, most of the polynomials in the result have the same monomials for their support, if p varies. This makes it very likely that the correct representative over ℚ has those same monomials as support and that the polynomials we computed mod p are the mod-p versions of the desired polynomials over ℚ.</p>
<p>This can fail for some primes which we need to ignore. In the function <code>likely_support</code>, we use a simple ad-hoc heuristic to find the exponents that probably occur in the ℚ-result: if a monomial occurs in 50% of the mod-p polynomials (i.e. modulo 50% of the primes or more), we expect that it should occur in the ℚ-result as well. The 50% threshold does not matter much. We could have picked 10% or 90% instead, I suppose.</p>
<p>We only keep the primes p for which all the polynomials in the lift of 1 mod p have the expected supports. Out of the 2800 primes we started with, 2358 remain.</p>
<p>Next, we use the Chinese remainder theorem <code>CRT_list</code> to lift the coefficients to ℤ modulo the product of the primes and then <code>rational_reconstruction</code> to obtain rational coefficients, both implemented in <code>lift_poly_crt</code>. For <code>rational_reconstruction</code> to work, the product of the primes must have about twice as many digits as the numerators and denominators in the result over ℚ (which I do not know how to determine in advance). In our case, we obtain a modulus with 10538 digits.</p>
<p>I have used <code>cached_function.precompute</code> in order to parallelize some computations because it is a bit less awkward to use than the <code>@parallel</code> decorator, but also because it allows to increase the number of primes until rational reconstruction works, without redoing all the previous slow computations.</p>
<p>The slowest part of the computation is lifting 1 modulo the 2800 primes, so this benefits a lot from parallelization. Each prime takes about 10 seconds, so this would take about 7.5 hours on a single thread.</p>
<pre><code>R1.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
I1=R1.ideal(A1)
C1=list(I1.groebner_basis())
R2.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
I2=R2.ideal(A2)
C2=list(I2.groebner_basis())
R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
I=R.ideal(C)
########################################################################
@cached_function
def lift_one_modp(Ip):
Fp = Ip.ring()(1).lift(Ip.ideal())
print("completed p = %s" % Fp.ring().characteristic())
return Fp
@cached_function
def lift_poly_crt(fps):
ps = [fp.parent().characteristic() for fp in fps]
modulus = prod(ps)
Ds = [fp.dict() for fp in fps]
supp = set([e for Di in Ds for e in Di])
return R({e: CRT_list([Di[e].lift() if e in Di else 0 for Di in Ds],
list(ps)).rational_reconstruction(modulus) for e in supp})
def likely_support(fps):
supps = [set(fp.dict().keys()) for fp in fps]
supp_union = set.union(*supps)
supp = set()
half = len(supps) * 1/2 # exact amount does not matter
for e in supp_union:
k = 0
for s in supps:
if e in s:
k += 1
if k > half: # if exponent vector is in more than half the supports, we will probably want to keep it
supp.add(e)
break
return supp
ps = sorted(RecursivelyEnumeratedSet([next_prime(16000)], lambda p: [next_prime(p)])[:2800]) # many primes
Rps = [R.change_ring(GF(p)) for p in ps]
Ips = [I.change_ring(Rp).gens() for Rp in Rps]
CPUs = 16 # choose based on hardware
lift_one_modp.precompute(Ips, num_processes=CPUs)
assert all(type(r) != str for r in lift_one_modp.cache.values()) # otherwise, there was an exception
Fps = [lift_one_modp(Ip) for Ip in Ips]
supps = [likely_support(fps) for fps in zip(*Fps)]
# keep indices of primes for which all polynomials have the expected support
idxs = [j for j, Fp in enumerate(Fps) if all(set(fp.dict()) == supp for fp, supp in zip(Fp, supps))]
# restrict the results to those primes
ps1 = tuple(ps[j] for j in idxs)
Fps1 = [Fps[j] for j in idxs]
fpss = [tuple(fps) for fps in zip(*Fps1)]
modulus1 = prod(ps1)
floor(log(modulus1, 10)) # 10538 (must be large enough for rational reconstruction: this needs about twice as many digits as the numerators or denominators in the result have)
lift_poly_crt.precompute([(fps,) for fps in fpss], num_processes=CPUs)
assert all(type(r) != str for r in lift_poly_crt.cache.values()) # otherwise, there was an exception (if so, remove it from cache by lift_poly_crt.cache.clear())
F = [lift_poly_crt(fps) for fps in fpss] # the final result</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57247#post-id-57247Very interesting and useful answer!Mon, 24 May 2021 04:25:06 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57247#post-id-57247Comment by rburing for <p>Note that computing the lift of 1 with respect to your <code>C</code> modulo a prime is fast – it takes only a few seconds. Doing this for many primes can allow to lift the result to ℚ. In this case, the result <code>F</code> (<s>see attachment</s> <a href="https://myshare.uni-osnabrueck.de/f/cfe7fa26c18d47c691c8/?dl=1">uploaded here (31MB)</a>, for now, as the file seems too large for an attachment) is huge, though, with coefficients as large as 4826 and 4818 digits for numerator and denominator.</p>
<pre><code>sage: F, C = load("lift_one_result_F_C.sobj") # saved with Sage 9.3
sage: sum(f*c for f,c in zip(F, C)) # 1 minute
1</code></pre>
<p>I am not sure if this solution would be unique. Maybe there is a smaller one.</p>
<p>To obtain the result with respect to <code>A1</code> and <code>A2</code>, you could try to lift <code>C1</code> and <code>C2</code> individually, but this might make the result even larger.</p>
<hr>
<h2>Update</h2>
<p>When computing the lift of 1 mod p, most of the polynomials in the result have the same monomials for their support, if p varies. This makes it very likely that the correct representative over ℚ has those same monomials as support and that the polynomials we computed mod p are the mod-p versions of the desired polynomials over ℚ.</p>
<p>This can fail for some primes which we need to ignore. In the function <code>likely_support</code>, we use a simple ad-hoc heuristic to find the exponents that probably occur in the ℚ-result: if a monomial occurs in 50% of the mod-p polynomials (i.e. modulo 50% of the primes or more), we expect that it should occur in the ℚ-result as well. The 50% threshold does not matter much. We could have picked 10% or 90% instead, I suppose.</p>
<p>We only keep the primes p for which all the polynomials in the lift of 1 mod p have the expected supports. Out of the 2800 primes we started with, 2358 remain.</p>
<p>Next, we use the Chinese remainder theorem <code>CRT_list</code> to lift the coefficients to ℤ modulo the product of the primes and then <code>rational_reconstruction</code> to obtain rational coefficients, both implemented in <code>lift_poly_crt</code>. For <code>rational_reconstruction</code> to work, the product of the primes must have about twice as many digits as the numerators and denominators in the result over ℚ (which I do not know how to determine in advance). In our case, we obtain a modulus with 10538 digits.</p>
<p>I have used <code>cached_function.precompute</code> in order to parallelize some computations because it is a bit less awkward to use than the <code>@parallel</code> decorator, but also because it allows to increase the number of primes until rational reconstruction works, without redoing all the previous slow computations.</p>
<p>The slowest part of the computation is lifting 1 modulo the 2800 primes, so this benefits a lot from parallelization. Each prime takes about 10 seconds, so this would take about 7.5 hours on a single thread.</p>
<pre><code>R1.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6>=PolynomialRing(QQ,10)
A1=[u0 + 7/5*u1 + 7/5*u2 - 4/125,
5*v0 + 5*v1 + 7*v3 + 7*v5 + 1/5,
25*v0^2 + 25*v1^2 + 35*v3^2 + 35*v5^2 - 4/5,
5*v0^3 + 5*v1^3 + 7*v3^3 + 7*v5^3 - v0^2 + 1/125,
5*v0*v1^2 + 5*v1*v2^2 + 7*v3*v4^2 + 7*v5*v6^2 + 1/125,
5*u0*v1 - v1^2 + 7*u1*v3 + 7*u2*v5 + 1/125,
5*v1 + 5*v2 + 7*v4 + 7*v6 + 1/5,
25*v0*v1 + 25*v1*v2 + 35*v3*v4 + 35*v5*v6 + 1/5,
5*v0^2*v1 + 5*v1^2*v2 + 7*v3^2*v4 + 7*v5^2*v6 - v1^2 + 1/125,
25*v1^2 + 25*v2^2 + 35*v4^2 + 35*v6^2 - 4/5,
5*v1^3 + 5*v2^3 + 7*v4^3 + 7*v6^3 - u0 + 1/125,
5*u0*v2 - v2^2 + 7*u1*v4 + 7*u2*v6 + 1/125]
I1=R1.ideal(A1)
C1=list(I1.groebner_basis())
R2.<y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,10)
A2=[y0 + 7/5*y1 + 7/5*y2 - 4/125,
5*z0 + 5*z1 + 7*z3 + 7*z5 + 1/5,
25*z0^2 + 25*z1^2 + 35*z3^2 + 35*z5^2 - 4/5,
5*z0^3 + 5*z1^3 + 7*z3^3 + 7*z5^3 - y0 + 1/125,
5*y0*z0 - z0^2 + 7*y1*z3 + 7*y2*z5 + 1/125,
5*z0*z1^2 + 5*z1*z2^2 + 7*z3*z4^2 + 7*z5*z6^2 - z1^2 + 1/125,
5*z1 + 5*z2 + 7*z4 + 7*z6 + 1/5,
25*z0*z1 + 25*z1*z2 + 35*z3*z4 + 35*z5*z6 + 1/5,
5*z0^2*z1 + 5*z1^2*z2 + 7*z3^2*z4 + 7*z5^2*z6 + 1/125,
5*y0*z1 - z1^2 + 7*y1*z4 + 7*y2*z6 + 1/125,
25*z1^2 + 25*z2^2 + 35*z4^2 + 35*z6^2 - 4/5,
5*z1^3 + 5*z2^3 + 7*z4^3 + 7*z6^3 - z2^2 + 1/125]
I2=R2.ideal(A2)
C2=list(I2.groebner_basis())
R.<u0,u1,u2,v0,v1,v2,v3,v4,v5,v6,y0,y1,y2,z0,z1,z2,z3,z4,z5,z6>=PolynomialRing(QQ,20)
C=C1+C2+[5*u0*z2 + 5*u1*z4 + 7*u2*z6 - u0 + 1/125] # the last one is the extra one
I=R.ideal(C)
########################################################################
@cached_function
def lift_one_modp(Ip):
Fp = Ip.ring()(1).lift(Ip.ideal())
print("completed p = %s" % Fp.ring().characteristic())
return Fp
@cached_function
def lift_poly_crt(fps):
ps = [fp.parent().characteristic() for fp in fps]
modulus = prod(ps)
Ds = [fp.dict() for fp in fps]
supp = set([e for Di in Ds for e in Di])
return R({e: CRT_list([Di[e].lift() if e in Di else 0 for Di in Ds],
list(ps)).rational_reconstruction(modulus) for e in supp})
def likely_support(fps):
supps = [set(fp.dict().keys()) for fp in fps]
supp_union = set.union(*supps)
supp = set()
half = len(supps) * 1/2 # exact amount does not matter
for e in supp_union:
k = 0
for s in supps:
if e in s:
k += 1
if k > half: # if exponent vector is in more than half the supports, we will probably want to keep it
supp.add(e)
break
return supp
ps = sorted(RecursivelyEnumeratedSet([next_prime(16000)], lambda p: [next_prime(p)])[:2800]) # many primes
Rps = [R.change_ring(GF(p)) for p in ps]
Ips = [I.change_ring(Rp).gens() for Rp in Rps]
CPUs = 16 # choose based on hardware
lift_one_modp.precompute(Ips, num_processes=CPUs)
assert all(type(r) != str for r in lift_one_modp.cache.values()) # otherwise, there was an exception
Fps = [lift_one_modp(Ip) for Ip in Ips]
supps = [likely_support(fps) for fps in zip(*Fps)]
# keep indices of primes for which all polynomials have the expected support
idxs = [j for j, Fp in enumerate(Fps) if all(set(fp.dict()) == supp for fp, supp in zip(Fp, supps))]
# restrict the results to those primes
ps1 = tuple(ps[j] for j in idxs)
Fps1 = [Fps[j] for j in idxs]
fpss = [tuple(fps) for fps in zip(*Fps1)]
modulus1 = prod(ps1)
floor(log(modulus1, 10)) # 10538 (must be large enough for rational reconstruction: this needs about twice as many digits as the numerators or denominators in the result have)
lift_poly_crt.precompute([(fps,) for fps in fpss], num_processes=CPUs)
assert all(type(r) != str for r in lift_poly_crt.cache.values()) # otherwise, there was an exception (if so, remove it from cache by lift_poly_crt.cache.clear())
F = [lift_poly_crt(fps) for fps in fpss] # the final result</code></pre>
https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57241#post-id-57241Nice work! The common denominator is divisible by the 785 digit number I found in the comments. The coefficients of the lift of 1 are not uniquely determined because you can always add syzygies (lifts of 0), and conversely if you have two different lifts of 1 then their difference is a syzygy.Sun, 23 May 2021 16:50:10 +0200https://ask.sagemath.org/question/57180/quick-trivial-groebner-basis-but-too-long-lift-of-one/?comment=57241#post-id-57241