gh-40810: Implement is_hyperelliptic
    
<!-- ^ Please provide a concise and informative title. -->
<!-- ^ Don't put issue numbers in the title, do this in the PR
description below. -->
<!-- ^ For example, instead of "Fixes #12345" use "Introduce new method
to calculate 1 + 2". -->
<!-- v Describe your changes below in detail. -->
<!-- v Why is this change required? What problem does it solve? -->
<!-- v If this PR resolves an open issue, please link to it here. For
example, "Fixes #12345". -->

This PR implements a method to determine if a Riemann surface is
hyperelliptic. Such a method [exists in Maple](https://www.maplesoft.com
/support/help/Maple/view.aspx?path=algcurves/is_hyperelliptic) and this
provides an open-source alternative following the original method of
Maple. It keep Sage more up to date with the tools required for Riemann
surface research.

### 📝 Checklist

<!-- Put an `x` in all the boxes that apply. -->

- [x] The title is concise and informative.
- [x] The description explains in detail what this PR is about.
- [ ] I have linked a relevant issue or discussion.
- [x] I have created tests covering the changes.
- [x] I have updated the documentation and checked the documentation
preview.

### ⌛ Dependencies

<!-- List all open PRs that this PR logically depends on. For example,
-->
<!-- - #12345: short description why this is a dependency -->
<!-- - #34567: ... -->
    
URL: #40810
Reported by: Linden Disney-Hogg
Reviewer(s): Linden Disney-Hogg, Michael Orlitzky

Author: Release Manager

src/doc/en/reference/references/index.rst

@@ -2424,6 +2424,13 @@ REFERENCES:
 
 .. [dotspec] http://www.graphviz.org/doc/info/lang.html
 
+.. [DP2011] \B. Deconinck and M. S. Patterson, *Computing with plane algebraic
+            curves and Riemann surfaces: The algorithms of the Maple package
+            "Algcurves"*, In: A. Bobenko and C. Klein (eds) Computational 
+            approach to Riemann surfaces. Lecture Notes in Mathematics 2013. 
+            Springer, Berlin, Heidelberg. (2011).
+            :doi:`10.1007/978-3-642-17413-1_2`
+
 .. [DPS2017] Kevin Dilks, Oliver Pechenik, and Jessica Striker,
              *Resonance in orbits of plane partitions and increasing
              tableaux*, JCTA 148 (2017), 244-274,

src/sage/schemes/riemann_surfaces/riemann_surface.py

@@ -1908,6 +1908,66 @@ def cohomology_basis(self, option=1):
             self._differentials = generators
         return self._differentials
 
+    def is_hyperelliptic(self):
+        r"""
+        Determine whether ``self`` is a hyperelliptic curve.
+
+        Uses the subspace of quadratic differentials generated by the
+        holomorphic differentials on ``self`` to determine if the curve is
+        hyperelliptic. Elliptic curves, i.e. curves of genus 1, are considered
+        to be hyperelliptic.
+
+        OUTPUT:
+
+        - boolean; whether or not the curve is hyperelliptic
+
+        EXAMPLES:
+
+        We consider the three examples given in [DP2011]_::
+
+            sage: from sage.schemes.riemann_surfaces.riemann_surface import RiemannSurface
+            sage: R.<x,y> = QQ[]
+            sage: f1 = y^3 - 2*x*y + x^4
+            sage: S1 = RiemannSurface(f1)
+            sage: S1.is_hyperelliptic()
+            True
+            sage: f2 = y^3 - 2*x*y + x^5
+            sage: S2 = RiemannSurface(f2)
+            sage: S2.is_hyperelliptic()
+            False
+            sage: f3 = y^9 + 3*x^2*y^6 + 3*x^4*y^3 + x^6 + y^2
+            sage: S3 = RiemannSurface(f3) # long time
+            sage: S3.is_hyperelliptic() # long time
+            True
+
+        ALGORITHM:
+
+        This uses the approach of van Hoeij who implemented the same function
+        in ``algcurves``. This uses the theory of Noether, and is described in
+        [DP2011]_. It finds the rank of the subspace of quadratic differentials
+        generated by holomorphic differentials: if this rank is `< 2g`, then
+        the curve is hyperelliptic.
+        """
+        # This uses common theory of Riemann surfaces
+        if self.genus <= 2:
+            return True
+
+        # Find the products of pairs of holomorphic differentials. We reduce
+        # these modulo the curve so that a linear dependence between the
+        # monomial coefficients becomes a linear dependence between the
+        # corresponding quadratic differentials
+        II = self._R.ideal(self.f)
+        diffs = self.cohomology_basis()
+        pairs = [II.reduce(diffs[i]*diffs[j]) for i in range(self.genus)
+                 for j in range(i, self.genus)]
+
+        # Find the monomials present and their coefficients
+        mons = {mon for p in pairs for mon in p.monomials()}
+        CM = Matrix([[p.monomial_coefficient(mon) for p in pairs]
+                     for mon in mons])
+        # test the number of linearly independent pairs
+        return CM.rank() <= 2*self.genus - 1
+
     def _bounding_data(self, differentials, exact=False):
         r"""
         Compute the data required to bound a differential on a circle.