gh-35991: Implement Goss polynomials of Drinfeld modules
    
<!-- ^^^^^
Please provide a concise, informative and self-explanatory title.
Don't put issue numbers in there, do this in the PR body below.
For example, instead of "Fixes #1234" use "Introduce new method to
calculate 1+1"
-->
<!-- Describe your changes here in detail -->

<!-- Why is this change required? What problem does it solve? -->
<!-- If this PR resolves an open issue, please link to it here. For
example "Fixes #12345". -->
<!-- If your change requires a documentation PR, please link it
appropriately. -->

This PR implements Goss polynomials of Drinfeld modules. The interested
PR reviewer should read section 3 of [this
paper](https://link.springer.com/article/10.1007/BF01410204). In
particular, the recurrence relation (ii) of Proposition 3.4 of op. cit.
is used in order to compute the $n$-th Goss polynomial of a Drinfeld
module.

CC: @xcaruso @ymusleh @kryzar

### 📝 Checklist

<!-- Put an `x` in all the boxes that apply. -->
<!-- If your change requires a documentation PR, please link it
appropriately -->
<!-- If you're unsure about any of these, don't hesitate to ask. We're
here to help! -->
<!-- Feel free to remove irrelevant items. -->

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

### ⌛ Dependencies

<!-- List all open PRs that this PR logically depends on
- #12345: short description why this is a dependency
- #34567: ...
-->

<!-- If you're unsure about any of these, don't hesitate to ask. We're
here to help! -->
    
URL: #35991
Reported by: David Ayotte
Reviewer(s): Antoine Leudière, Xavier Caruso

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

@@ -2665,6 +2665,9 @@ REFERENCES:
             TR-737-05, (2005).
             ftp://ftp.cs.princeton.edu/reports/2005/737.pdf
 
+.. [Gek1988] \E.-U. Gekeler, On the coefficients of Drinfel'd modular
+             forms. Invent. Math. 93 (1988), no. 3, 667-700
+
 .. [Gek1991] \E.-U. Gekeler. On finite Drinfeld modules. Journal of
              algebra, 1(141):187–203, 1991.
 

src/sage/rings/function_field/drinfeld_modules/drinfeld_module.py

@@ -1361,6 +1361,102 @@ def gen(self):
         """
         return self._gen
 
+    @cached_method
+    def _compute_goss_polynomial(self, n, q, poly_ring, X):
+        r"""
+        Utility function for computing the n-th Goss polynomial.
+
+        The user should not call this method directly, but
+        :meth:`goss_polynomial` instead.
+
+        TESTS::
+
+            sage: A = GF(2^2)['T']
+            sage: K.<T> = Frac(A)
+            sage: phi = DrinfeldModule(A, [T, T+1, T^2, 1])
+            sage: poly_ring = phi.base()['X']
+            sage: X = poly_ring.gen()
+            sage: phi._compute_goss_polynomial(0, 2^2, poly_ring, X)
+            0
+            sage: phi._compute_goss_polynomial(3, 2^2, poly_ring, X)
+            X^3
+            sage: phi._compute_goss_polynomial(4*3, 2^2, poly_ring, X)
+            X^12
+            sage: phi._compute_goss_polynomial(9, 2^2, poly_ring, X)
+            X^9 + (1/(T^3 + T^2 + T))*X^6 + (1/(T^6 + T^4 + T^2))*X^3
+
+        """
+        # Trivial cases
+        if n.is_zero():
+            return poly_ring.zero()
+        if n <= q - 1:
+            return X**n
+        if n % q == 0:
+            return self.goss_polynomial(n // q)**q
+        # General case
+        pol = sum(self._compute_coefficient_exp(i+1)
+                  *self._compute_goss_polynomial(n - q**(i+1), q, poly_ring, X)
+                  for i in range(0, (n.log(q).n()).floor()))
+        return X*(self._compute_goss_polynomial(n - 1, q, poly_ring, X) + pol)
+
+    def goss_polynomial(self, n, var='X'):
+        r"""
+        Return the `n`-th Goss polynomial of the Drinfeld module.
+
+        Note that Goss polynomials are only defined for Drinfeld modules
+        of characteristic zero.
+
+        INPUT:
+
+        - ``n`` (integer) -- the index of the Goss polynomial
+
+        - ``var`` (str, default: ``'X'``) -- the name of polynomial
+          variable.
+
+        OUTPUT:
+
+        - a univariate polynomial in ``var`` over the base `A`-field.
+
+        EXAMPLES::
+
+            sage: A = GF(3)['T']
+            sage: K.<T> = Frac(A)
+            sage: phi = DrinfeldModule(A, [T, 1])  # Carlitz module
+            sage: phi.goss_polynomial(1)
+            X
+            sage: phi.goss_polynomial(2)
+            X^2
+            sage: phi.goss_polynomial(4)
+            X^4 + (1/(T^3 + 2*T))*X^2
+            sage: phi.goss_polynomial(5)
+            X^5 + (2/(T^3 + 2*T))*X^3
+            sage: phi.goss_polynomial(10)
+            X^10 + (1/(T^3 + 2*T))*X^8 + (1/(T^6 + T^4 + T^2))*X^6 + (1/(T^9 + 2*T^3))*X^4 + (1/(T^18 + 2*T^12 + 2*T^10 + T^4))*X^2
+
+        TESTS::
+
+            sage: Fq.<z> = GF(25)
+            sage: A.<T> = Fq[]
+            sage: phi = DrinfeldModule(A, [z, 1])
+            sage: phi.goss_polynomial(1)
+            Traceback (most recent call last):
+            ...
+            ValueError: characteristic must be zero (=T^2 + 4*T + 2)
+
+        REFERENCE:
+
+        Section 3 of [Gek1988]_ provides an exposition of Goss
+        polynomials.
+        """
+        if self.category()._characteristic:
+            raise ValueError(f"characteristic must be zero (={self.characteristic()})")
+        n = ZZ(n)
+        K = self.base()
+        poly_ring = K[var]
+        X = poly_ring.gen()
+        q = self._Fq.cardinality()
+        return self._compute_goss_polynomial(n, q, poly_ring, X)
+
     def height(self):
         r"""
         Return the height of the Drinfeld module if the function field