sage: Fix incompatibility with sage 9.4

real-or-random commented at 11:16 am on January 31, 2022: contributor

allexprs is already the product all numerators. Don’t take it’s numerator again.

To be honest, I have no idea what I’m doing here.

`allexprs` is already the product all numerators. Don't take it's
numerator again.

Fixes #1067.

e108d0039c

jonasnick commented at 1:43 pm on January 31, 2022: contributor

I get the same errors here as with the change_ring “workaround”:

0Formula secp256k1_gej_add_var:
1  add:
2    branch 4: FAIL, <map object at 0x7f417d48b820> fails (assuming )
3
4Formula secp256k1_gej_add_zinv_var:
5  add:
6    branch 4: FAIL, <map object at 0x7f417d67d420> fails (assuming )

sage: Fix printing of errors

Python 3 often returns iterable map objects where Python 2 returned
list. We can just them down to lists explicitly.

Overlooked in 13c88efed0005eb6745a222963ee74564054eafb.

b54d843eac

real-or-random commented at 5:00 pm on January 31, 2022: contributor

@jonasnick Oh I missed this… I believe this a real of case of incompleteness (something in sage changed?) and now our prover isn’t complete enough to prove the formula. But the fix is easy. If you prints() to the relevant part of prove_nonzero like this…

0  print("nonzero", nonzero)
1  for (f, n) in allexprs.factor():
2    print("factor", f)

… you get

0nonzero {by, ax, bx, Az, bx - ax, ay, Bz}
1factor -bx + ax

So it put bx - ax in the set of nonzero expressions but then fails to see that then -bx + ax is nonzero as well. I added a commit to address this.

Disclaimer: I don’t really know what this code is doing in the end, so please take this with a grain of salt. @sipa should check if I’m right.

sipa cross-referenced this on Jan 31, 2022 from issue TypeError in prove_group_implementations.sage by jonasnick

in sage/group_prover.sage:225 in eb86372d5d outdated

221@@ -222,12 +222,12 @@ def prove_nonzero(R, exprs, assume):
222       return (False, [exprs[expr]])
223   allexprs = reduce(lambda a,b: numerator(a)*numerator(b), exprs, 1)
224   for (f, n) in allexprs.factor():
225-    if f not in nonzero:
226+    if f not in nonzero and -f not in nonzero:

sipa commented at 5:21 pm on January 31, 2022:

I wonder if this is generic enough, and there couldn’t be arbitrary constant factor differences. If the results aren’t normalized in a way that result in unique signs, other factors might differ too. What about:

0  allexprs = reduce(lambda a,b: numerator(a)*numerator(b), exprs, 1)
1  for (f, n) in allexprs.factor():
2    if not any((f % nz) == 0 and (nz % f) == 0 for nz in nonzero):
3      ok = False

real-or-random commented at 6:40 pm on January 31, 2022:

Both the elements of nonzero and of allexprs.factor() are prime factors obtained by factorization, so I think they can differ at most by their unit() (https://doc.sagemath.org/html/en/reference/structure/sage/structure/factorization.html#sage.structure.factorization.Factorization.unit), which is -1 or 1.

So an element of nonzero or allexprs.factor() cannot have an additional factor c with c != -1 and c != 1 because then it wouldn’t be a prime factor. But there could be a constant prime element c in allexprs.factor(), which is then trivially nonzero but the prover fails this.

So if we want to allow for this, I’d suggest (which does not on the % operator for which I can’t find docs on):

0  for (f, n) in allexprs.factor():
1    if f not in nonzero and -f not in nonzero and not f.is_constant():
2      ok = False

Docs are here: https://doc.sagemath.org/html/en/reference/polynomial_rings/sage/rings/polynomial/multi_polynomial_element.html#sage.rings.polynomial.multi_polynomial_element.MPolynomial_polydict.is_constant

sipa commented at 6:43 pm on January 31, 2022:

Are you sure? These are (multivariate) polynomial rings, and all constant terms are units for that.

Specifically, I’m concerned about the case where you’d have factors a + 2*b or b + 1/2*a, depending on which variable is considered the one to make monic.

real-or-random commented at 7:23 pm on January 31, 2022:

Are you sure? These are (multivariate) polynomial rings, and all constant terms are units for that.

I’m not sure at all.

Specifically, I’m concerned about the case where you’d have factors a + 2b or b + 1/2a, depending on which variable is considered the one to make monic.

Hm yeah, I’d imagine this is the same everytime in the same ring class in sage but I don’t know.

Would you then change every similar if within this function?

sipa commented at 7:27 pm on January 31, 2022:

Hm yeah, I’d imagine this is the same everytime in the same ring class in sage but I don’t know.

Indeed, but the same reasoning should IMO apply to the just the factor 1 or -1. ~~If there is a canonical ordering of factors, and e.g. factor() always returns the lexicographically smallest list of factors, there should never be any ambiguity at all.~~ The fact that the sign differs makes me suspicious, so I’d rather be overly cautious (the code still works fine with that change).

EDIT: I realize the above is incorrect. Just because individual inputs to factor() are canonical doesn’t imply the same is true for all factors individually. Still, I’m unsure if 1 and -1 are the only possible differences.

Would you then change every similar if within this function?

Yeah. Perhaps introduce a unit_ratio_between(a, b) function or so to abstract it and make it a bit more readable?

real-or-random commented at 11:46 am on February 1, 2022:

Ok I played around and I’m 99% convinced that 1 and -1 are complete. The reason is that the unit part is ignored when we convert the Factorization to a list in our code (one of the first examples in https://doc.sagemath.org/html/en/reference/structure/sage/structure/factorization.html), which then gives a canonical representation.

Demonstration: (Sorry now even sage 9.5… Arch Linux is pretty quick with its updates…)

 0┌────────────────────────────────────────────────────────────────────┐
 1│ SageMath version 9.5, Release Date: 2022-01-30                     │
 2│ Using Python 3.10.2. Type "help()" for help.                       │
 3└────────────────────────────────────────────────────────────────────┘
 4sage: R.<ax,bx,ay,by,Az,Bz,Ai,Bi> = PolynomialRing(QQ,8,order='invlex')
 5sage: R
 6Multivariate Polynomial Ring in ax, bx, ay, by, Az, Bz, Ai, Bi over Rational Field
 7
 8sage: R(10 * ax + bx / 2).factor()
 9(1/2) * (bx + 20*ax)
10sage: R(20 * ax + bx).factor()
11bx + 20*ax
12
13sage: R(20 * ax + bx).factor().unit()
141
15sage: R(10 * ax + bx/2).factor().unit()
161/2
17
18sage: list(R(20 * ax + bx).factor())
19[(bx + 20*ax, 1)]
20sage: list(R(10 * ax + bx / 2).factor())
21[(bx + 20*ax, 1)]

Why do 1 and -1 matter then at all, you may ask?

0sage: R((-2 * (bx - ax)) ^ 1).factor()
1(-2) * (bx - ax)
2sage: R((-2 * (bx - ax)) ^ 2).factor()
3(4) * (-bx + ax)^2
4sage: R((-2 * (bx - ax)) ^ 3).factor()
5(8) * (-bx + ax)^3

It’s not very consistent where they put the minus: it’s not even consistent among odd exponents.

And indeed, I debugged the code on master and we hit (bx - ax)^34 in allexprs, where the factor then gets a sign change during factorization.

For curiosity, can you check what this code does on your Sage 9.2?

The canonical representation (modulo the sign) is such that all fractional parts are moved to the unit:

0sage: R((1/10 * ax + 1/15*bx)*(ax/7+bx)).factor()
1(1/210) * (2*bx + 3*ax) * (7*bx + ax)

Of course, this might change in the future but all of sage may change in the future, so I think we should provide a fix only when it actually breaks.

So with all that said, I’m in favor of sticking to the current solution in this PR.

jonasnick commented at 10:19 pm on February 2, 2022:

I did some independent experiments and the only factorizations I was able to get Sage (9.4) to produce had the coefficients of the non-unit factors be integers with gcd 1. This would make them unique if it weren’t for this +,- inconsistency.

I considered that this is a regression in Sage >9.2 with order="invlex", because replacing it with order="lex" makes the script work without the -f not in nonzero fix. But no.

0sage: a, b, c = QQ["a,b,c"].gens()
1....: f = 1 / 2 * a - b
2....: f.factor()
3....: 
4(1/2) * (a - 2*b)
5sage: f2 = (1 / 2 * a - b) * (a - c)
6....: f2.factor()
7(1/2) * (-a + c) * (-a + 2*b)

I’d be fine with leaving the PR as is. Another solution would be to wrap the factor() method such that it behaves as we would expect, i.e. the smallest variable name in a factor gets a positive coefficient (and our other assumptions are asserted).

jonasnick commented at 10:23 pm on February 2, 2022:

nit: perhaps note in the commit description that this commit is required to make the script pass with sage 9.4.

sipa commented at 10:27 pm on February 2, 2022:

 0$ sage
 1┌────────────────────────────────────────────────────────────────────┐
 2│ SageMath version 9.2, Release Date: 2020-10-24                     │
 3│ Using Python 3.9.7. Type "help()" for help.                        │
 4└────────────────────────────────────────────────────────────────────┘
 5sage: R.<ax,bx,ay,by,Az,Bz,Ai,Bi> = PolynomialRing(QQ,8,order='invlex')
 6sage: R
 7Multivariate Polynomial Ring in ax, bx, ay, by, Az, Bz, Ai, Bi over Rational Field
 8sage: R(10 * ax + bx / 2).factor()
 9(1/2) * (bx + 20*ax)
10sage: R(20 * ax + bx).factor()
11bx + 20*ax
12sage: R(20 * ax + bx).factor().unit()
131
14sage: R(10 * ax + bx/2).factor().unit()
151/2
16sage: list(R(20 * ax + bx).factor())
17[(bx + 20*ax, 1)]
18sage: list(R(10 * ax + bx / 2).factor())
19[(bx + 20*ax, 1)]
20sage: R((-2 * (bx - ax)) ^ 1).factor()
21(2) * (-bx + ax)
22sage: R((-2 * (bx - ax)) ^ 2).factor()
23(4) * (-bx + ax)^2
24sage: R((-2 * (bx - ax)) ^ 3).factor()
25(8) * (-bx + ax)^3
26sage: R((1/10 * ax + 1/15*bx)*(ax/7+bx)).factor()
27(1/210) * (2*bx + 3*ax) * (7*bx + ax)

sipa commented at 10:41 pm on February 2, 2022:

Is this (f not in nonzero) and (-f not in nonzero) needed in other places where this pattern occurs below?

real-or-random commented at 12:02 pm on February 3, 2022:

Another solution would be to wrap the factor() method such that it behaves as we would expect, i.e. the smallest variable name in a factor gets a positive coefficient (and our other assumptions are asserted).

That’s a clean idea. I did this.

Is this (f not in nonzero) and (-f not in nonzero) needed in other places where this pattern occurs below?

Not for the current group laws but I added the wrapper to the other places where we call factor().

real-or-random force-pushed on Feb 3, 2022

real-or-random commented at 12:02 pm on February 3, 2022: contributor

I added sage to CI. Ready from my side (if CI passes…)

real-or-random commented at 1:46 pm on February 3, 2022: contributor

Upstream report: https://trac.sagemath.org/ticket/33283

in sage/group_prover.sage:203 in bed969b04f outdated

198+  # Assert p is not 0 and that its non-zero coeffients are coprime.
199+  # (We could just work with the primitive part p/p.content() but we want to be
200+  # aware if factor() does not return a primitive part in future sage versions.)
201+  assert p.content() == 1
202+  # Ensure that the first non-zero coefficient is positive.
203+  return p if p.lc() > 0 else -p

sipa commented at 7:45 pm on February 3, 2022:

What about using p / p.lc()? That will always guarantee a result with leading coefficient 1, even if there is more than 1/-1 variability.

jonasnick commented at 8:09 pm on February 3, 2022:

That seems possible but by asserting that the content is 1 we shouldn’t have more than 1/-1 variability.

sipa commented at 9:21 pm on February 3, 2022:

Sure, but if you use p / p.lc(), there is no need for such an assertion.

jonasnick commented at 9:58 pm on February 3, 2022:

Ok. It seemed slightly more appealing to me to just patch the sage bug but *shrug*

real-or-random commented at 8:15 am on February 4, 2022:

What about using p / p.lc()? That will always guarantee a result with leading coefficient 1, even if there is more than 1/-1 variability.

That seems to work, though this changes the representation further. p/p.content() if p.lc() > 0 else -p/p.content() will be closer to what we have now (and return a primitive polynomial).

But yeah, see the comment. It was an intentional decision to add the assertion. I think if that assertion fails, it may be the case that sage did larger changes and we may want to investigate this then.

But I can change it if you prefer.

sipa commented at 3:17 pm on February 4, 2022:

No problem.

jonasnick commented at 8:34 am on February 4, 2022: contributor

ACK bed969b04fa292346edf5766e681a25d68b7a9d7

in sage/prove_group_implementations.sage:303 in bed969b04f outdated

308-    check_exhaustive_jacobian_weierstrass("secp256k1_gej_add_var", 0, 7, 5, formula_secp256k1_gej_add_var, 43)
309-    check_exhaustive_jacobian_weierstrass("secp256k1_gej_add_ge_var", 0, 7, 5, formula_secp256k1_gej_add_ge_var, 43)
310-    check_exhaustive_jacobian_weierstrass("secp256k1_gej_add_zinv_var", 0, 7, 5, formula_secp256k1_gej_add_zinv_var, 43)
311-    check_exhaustive_jacobian_weierstrass("secp256k1_gej_add_ge", 0, 7, 16, formula_secp256k1_gej_add_ge, 43)
312-    check_exhaustive_jacobian_weierstrass("secp256k1_gej_add_ge_old [should fail]", 0, 7, 4, formula_secp256k1_gej_add_ge_old, 43)
313+    success = success & check_exhaustive_jacobian_weierstrass("secp256k1_gej_add_var", 0, 7, 5, formula_secp256k1_gej_add_var, 43)

real-or-random commented at 2:27 pm on February 4, 2022:

This is buggy for --exhaustive, let me fix this.

real-or-random commented at 2:41 pm on February 4, 2022:

Ok done. Also cherry-picked a commit from https://github.com/real-or-random/secp256k1/commits/202202-gej_add_ge that I discovered when working on #1033.

jonasnick commented at 9:47 pm on February 4, 2022:

The commit looks fine but is there a specific reason why it was added?

real-or-random commented at 7:44 am on February 5, 2022:

I added it because it fixes a bug in the group prover that prevents us to look into other group formulas. In particular, it stopped me from playing around with the formula here because it may return a literal 1: https://github.com/real-or-random/secp256k1/blob/dc53e24a1b63ec07a9ab086280cab1355cd3a881/sage/prove_group_implementations.sage#L209

jonasnick commented at 9:57 pm on February 5, 2022:

Ah, thanks that makes sense.

sage: Exit with non-zero status in case of failures eae75869cf

sage: Normalize sign of polynomial factors in prover

The prover, when run on recent sage versions,  failed to prove some of its
goals due to a change in sage. This commit adapts our code accordingly.
The prover passes again after this commit.

77cfa98dbc

ci: Run sage prover on CI d8d54859ed

real-or-random force-pushed on Feb 4, 2022

sage: Ensure that constraints are always fastfracs

Even if they are constants created in the formula functions. We now
lift integer constants to fastfracs.

ebb1beea78

sipa commented at 4:21 pm on February 4, 2022: contributor

ACK ebb1beea7832207e5c8c5112d250fd216259ef41

jonasnick commented at 9:59 pm on February 5, 2022: contributor

ACK ebb1beea7832207e5c8c5112d250fd216259ef41

jonasnick merged this on Feb 5, 2022

jonasnick closed this on Feb 5, 2022

fanquake referenced this in commit 8f8631d826 on Mar 17, 2022

fanquake referenced this in commit 4bb1d7e62a on Mar 17, 2022

fanquake referenced this in commit 465d05253a on Mar 30, 2022

jonasnick cross-referenced this on Mar 30, 2022 from issue Upstream PRs 1064, 1049, 899, 1068, 1072, 1069, 1074, 1026, 1033, 748, 1079, 1088, 1090, 731, 1089, 995, 1094, 1093 by jonasnick

real-or-random referenced this in commit 6c0aecf72b on Apr 1, 2022

fanquake referenced this in commit afb7a6fe06 on Apr 6, 2022

gwillen referenced this in commit 35d6112a72 on May 25, 2022

patricklodder referenced this in commit 21badcf9d2 on Jul 25, 2022

patricklodder referenced this in commit 03002a9013 on Jul 28, 2022

janus referenced this in commit 3a0652a777 on Aug 4, 2022

real-or-random cross-referenced this on Feb 6, 2023 from issue Improve sign consistency in rational polynomial factorization by real-or-random

str4d referenced this in commit 522190d5c3 on Apr 21, 2023

vmta referenced this in commit e1120c94a1 on Jun 4, 2023

vmta referenced this in commit 8f03457eed on Jul 1, 2023

sage: Fix incompatibility with sage 9.4 #1068