Arithmetic on Jacobians of Hyperelliptic curves fails when there are two points at infinity

Description

This is a meta issue made from three encountered during Sage Days 123. It was noticed that there was ambiguous representation and failing group laws in certain scenarios. Issue #37093 noticed ambiguous representation and #37101 highlighted that the even degree model could not support the case for the even degree model due to a lack of data in the elements of the Jacobian.

Expected Behaviour

The proper fix for this is to rewrite the mumford coordinates and the arithmetic to allow for proper handling when there are two points at infinity. The Magma implementation of this follows: https://www.math.auckland.ac.nz/~sgal018/crypto-book/ch10.pdf and I think this is what SageMath should do too.

Issues

Group law fails when the model is even

For the even degree model, there are two points at infinity and the current implementation has no way of distinguishing between $P - (\infty+)$ and $P - (\infty-)$.

A proper fix for this would require refactoring how src/sage/schemes/hyperelliptic_curves/jacobian_morphism.py currently behaves for the even degree case. For instance some tuple of the form (u, v, ±). Magma handles this with a tuple (a, b, d).

I don't currently have a good idea on how to fix this, but I'll try and come up with a reasonable fix to propose. At the moment this issue is really a place holder.

Thanks to @sabrinakunzweiler for bringing this to my attention.

At the moment there is no way to represent both $P - (\infty+)$ and $P - (\infty-)$. This means that if we randomly sample divisors on the Jacobian we will have H.count_points() elements missing.

sage: F = GF(7)
sage: R.<x> = PolynomialRing(F)
sage: f = x^6 + x + 1
sage: H = HyperellipticCurve(f)
sage: J = H.jacobian()
sage: H.count_points()
[9]
sage: H.zeta_function().numerator()(1) # order of the Jacobian
67
sage: len(set(tuple(J.random_element()) for _ in range(2000)))
58
sage: assert 58 + 9 == 67

Group law fails when degree of h is too large

sage: F = GF(11)
....: R.<x> = PolynomialRing(F)
....: f = 6*x^5 + 6*x^4 + 9*x^3 + 9*x + 8
....: h = 5*x^3 + 3*x^2 + 5*x + 4
....: H = HyperellipticCurve(f, h)
....: J = H.jacobian()
....: D = J(H.lift_x(1))
....: o = H.zeta_function().numerator()(1) # order of Jacobian
....: o * D == J(0)
....: o * D
False
(x^2 + 4*x + 4, y + 4*x + 3)

Mumford coordinates are not uniquely represented

When summing divisors for even degree models, the cantor reduction sometimes does not finish the reduction process. This seems to be due to a TODO left in the code.

sage: F = GF(163)
sage: R.<x> = PolynomialRing(F)
sage: f = x^6 + x^2 + 1
sage: H = HyperellipticCurve(f)
sage: J = H.jacobian()
sage: Ps = H.points()
sage: Ds = [J(P) for P in Ps]
sage: sum(choice(Ds) for _ in range(7))
Returning ambiguous form of degree genus+1.
Returning ambiguous form of degree genus+1.
Returning ambiguous form of degree genus+1.
(x^3 + 11*x^2 + 104*x + 136, y + 95*x^2 + 144*x + 156)

Environment

- **OS**: MacOS
- **Sage Version**: 10.2

Checklist

[X] I have searched the existing issues for a bug report that matches the one I want to file, without success.
[X] I have read the documentation and troubleshoot guide

sagemath / sage