Fix PARI DCP solution.

3 files changed, 51 insertions, 26 deletions

diff --git a/pyecsca/sca/re/tree.py b/pyecsca/sca/re/tree.py
index 0f6bcf1..5d93f1d 100644
--- a/pyecsca/sca/re/tree.py
+++ b/pyecsca/sca/re/tree.py
@@ -212,6 +212,7 @@ def _build_tree(cfgs: Set[Any], *maps: Map, response: Optional[Any] = None) -> N
                 best_score = score
                 best_restricted = restricted
             # Early abort if optimal score is hit. The +1 is for "None" values which are not in the codomain.
+            # TODO: Move the None to the codomain and ditch the +1 as some codomains may not have it (complete).
             if score == ceil(n_cfgs / (len(dmap.codomain) + 1)):
                 break
     # We found nothing distinguishing the configs, so return them all (base case 2).
diff --git a/pyecsca/sca/re/zvp.py b/pyecsca/sca/re/zvp.py
index 98fcfa7..527c6ea 100644
--- a/pyecsca/sca/re/zvp.py
+++ b/pyecsca/sca/re/zvp.py
@@ -511,6 +511,11 @@ def solve_easy_dcp(xonly_polynomial: Poly, curve: EllipticCurve) -> Set[Point]:
     if has_pari:
         pari = cypari2.Pari()
         polynomial = pari(str(final.expr).replace("**", "^"))
+        if final.degree() == 0:
+            if final.is_zero:
+                # TODO: When the polynomial is zero, sympy does not give any roots, we want to return the generator.
+                pass
+            return set()
         roots = list(map(int, pari.polrootsmod(polynomial, curve.prime)))
     else:
         roots = final.ground_roots().keys()
@@ -531,40 +536,51 @@ def solve_hard_dcp(xonly_polynomial: Poly, curve: EllipticCurve, k: int) -> Set[
         # Find the roots (X1)
         roots = final.ground_roots().keys()
         # Finally lift the roots to find the points (if any)
+        # TODO: When the polynomial is zero, sympy does not give any roots, we want to return the generator.
     for root in roots:
         points.update(curve.affine_lift_x(Mod(int(root), curve.prime)))
     return points
 
 
 def solve_hard_dcp_cypari(xonly_polynomial: Poly, curve: EllipticCurve, k: int) -> Set[int]:
-    a, b = int(curve.parameters["a"]), int(curve.parameters["b"])
-    xonly_polynomial = subs_curve_params(xonly_polynomial, curve)
+    try:
+        a, b = int(curve.parameters["a"]), int(curve.parameters["b"])
+        xonly_polynomial = subs_curve_params(xonly_polynomial, curve)
 
-    # k^2 * degree
-    # k=25, deg=6, 128bit -> 3765, a 20MB
-    # k=32, deg=6, 128bit -> 6150, a 32MB
-    # k=10, deg=6, 128bit -> 606, a 4MB
-    outdegree = k**2 * xonly_polynomial.degree()
-    stacksize = 2 * (outdegree * (40 * curve.prime.bit_length()))
-    stacksizemax = 2 * stacksize
+        # k^2 * degree
+        # k=25, deg=6, 128bit -> 3765, a 20MB
+        # k=32, deg=6, 128bit -> 6150, a 32MB
+        # k=10, deg=6, 128bit -> 606, a 4MB
+        outdegree = k**2 * xonly_polynomial.degree()
+        # TODO: When can outdegree be 0?
+        # Magic heuristic, plus some constant term for very small polys
+        stacksize = 2 * (outdegree * (40 * curve.prime.bit_length())) + 1000000
+        stacksizemax = 15 * stacksize
 
-    pari = cypari2.Pari()
-    pari.allocatemem(stacksize, stacksizemax)
-    e = pari.ellinit([a, b], curve.prime)
-    mul = pari.ellxn(e, k)
-    x1, x2 = pari("x1"), pari("x2")
-    polynomial = pari(str(xonly_polynomial.expr).replace("**", "^"))
+        pari = cypari2.Pari()
+        pari.allocatemem(stacksize, stacksizemax, silent=True)
+        e = pari.ellinit([a, b], curve.prime)
+        mul = pari.ellxn(e, k)
+        x1, x2 = pari("x1"), pari("x2")
+        polynomial = pari(str(xonly_polynomial.expr).replace("**", "^"))
 
-    polydeg = pari.poldegree(polynomial, x2)
-    subspoly = 0
-    x = pari("x")
-    num, den = pari.subst(mul[0], x, x1), pari.subst(mul[1], x, x1)
-    for deg in range(polydeg+1):
-        monomial = pari.polcoef(polynomial, deg, x2)
-        monomial *= polypower(pari, num, deg)
-        monomial *= polypower(pari, den, polydeg-deg)
-        subspoly += monomial
-    res = set(map(int, pari.polrootsmod(subspoly, curve.prime)))
+        polydeg = pari.poldegree(polynomial, x2)
+        subspoly = 0
+        x = pari("x")
+        num, den = pari.subst(mul[0], x, x1), pari.subst(mul[1], x, x1)
+        for deg in range(polydeg+1):
+            monomial = pari.polcoef(polynomial, deg, x2)
+            monomial *= polypower(pari, num, deg)
+            monomial *= polypower(pari, den, polydeg-deg)
+            subspoly += monomial
+        if subspoly == pari.zero():
+            # TODO: Return the generator here when the API changes.
+            return set()
+        res = set(map(int, pari.polrootsmod(subspoly, curve.prime)))
+    except cypari2.PariError as err:
+        raise ValueError("PariError " + err.errtext())
+    except Exception as err:
+        raise ValueError(str(err))
     return res
 
 
diff --git a/test/sca/test_zvp.py b/test/sca/test_zvp.py
index 601a695..c590521 100644
--- a/test/sca/test_zvp.py
+++ b/test/sca/test_zvp.py
@@ -18,7 +18,7 @@ from pyecsca.sca.re.zvp import (
     zvp_points,
     compute_factor_set,
     addition_chain,
-    precomp_zvp_points,
+    precomp_zvp_points, solve_easy_dcp,
 )
 from pyecsca.ec.context import local, DefaultContext
 from sympy import symbols, Poly, sympify, FF
@@ -300,3 +300,11 @@ def test_big_boy(secp128r1, k):
     poly = Poly(poly_expr, domain=FF(secp128r1.curve.prime))
     res = zvp_points(poly, secp128r1.curve, k, secp128r1.order)
     assert res is not None
+
+
+@pytest.mark.parametrize("numero", [1, 0])
+def test_small_boy(secp128r1, numero):
+    x = symbols("x1")
+    poly = Poly(numero, x, domain=FF(secp128r1.curve.prime))
+    res = solve_easy_dcp(poly, secp128r1.curve)
+    assert res is not None
\ No newline at end of file