Drinfeld modules over rings of characteristic zero¶
This module provides the classes
sage.rings.function_fields.drinfeld_module.charzero_drinfeld_module.DrinfeldModule_charzero and
sage.rings.function_fields.drinfeld_module.charzero_drinfeld_module.DrinfeldModule_rational,
which both inherit
sage.rings.function_fields.drinfeld_module.drinfeld_module.DrinfeldModule.
AUTHORS:
- David Ayotte (2023-09) 
- Xavier Caruso (2024-12) - computation of class polynomials 
- class sage.rings.function_field.drinfeld_modules.charzero_drinfeld_module.DrinfeldModule_charzero(gen, category)[source]¶
- Bases: - DrinfeldModule- This class implements Drinfeld \(\mathbb{F}_q[T]\)-modules defined over fields of \(\mathbb{F}_q[T]\)-characteristic zero. - Recall that the \(\mathbb{F}_q[T]\)-characteristic is defined as the kernel of the underlying structure morphism. For general definitions and help on Drinfeld modules, see class - sage.rings.function_fields.drinfeld_module.drinfeld_module.DrinfeldModule.- Construction: - The user does not ever need to directly call - DrinfeldModule_charzero— the metaclass- DrinfeldModuleis responsible for instantiating the right class depending on the input:- sage: A = GF(3)['T'] sage: K.<T> = Frac(A) sage: phi = DrinfeldModule(A, [T, 1]) sage: phi Drinfeld module defined by T |--> t + T - >>> from sage.all import * >>> A = GF(Integer(3))['T'] >>> K = Frac(A, names=('T',)); (T,) = K._first_ngens(1) >>> phi = DrinfeldModule(A, [T, Integer(1)]) >>> phi Drinfeld module defined by T |--> t + T - sage: isinstance(phi, DrinfeldModule) True sage: from sage.rings.function_field.drinfeld_modules.charzero_drinfeld_module import DrinfeldModule_charzero sage: isinstance(phi, DrinfeldModule_charzero) True - >>> from sage.all import * >>> isinstance(phi, DrinfeldModule) True >>> from sage.rings.function_field.drinfeld_modules.charzero_drinfeld_module import DrinfeldModule_charzero >>> isinstance(phi, DrinfeldModule_charzero) True - Logarithm and exponential - It is possible to calculate the logarithm and the exponential of any Drinfeld modules of characteristic zero: - sage: A = GF(2)['T'] sage: K.<T> = Frac(A) sage: phi = DrinfeldModule(A, [T, 1]) sage: phi.exponential() z + ((1/(T^2+T))*z^2) + ((1/(T^8+T^6+T^5+T^3))*z^4) + O(z^8) sage: phi.logarithm() z + ((1/(T^2+T))*z^2) + ((1/(T^6+T^5+T^3+T^2))*z^4) + O(z^8) - >>> from sage.all import * >>> A = GF(Integer(2))['T'] >>> K = Frac(A, names=('T',)); (T,) = K._first_ngens(1) >>> phi = DrinfeldModule(A, [T, Integer(1)]) >>> phi.exponential() z + ((1/(T^2+T))*z^2) + ((1/(T^8+T^6+T^5+T^3))*z^4) + O(z^8) >>> phi.logarithm() z + ((1/(T^2+T))*z^2) + ((1/(T^6+T^5+T^3+T^2))*z^4) + O(z^8) - Goss polynomials - Goss polynomials are a sequence of polynomials related with the analytic theory of Drinfeld module. They provide a function field analogue of certain classical trigonometric functions: - sage: A = GF(2)['T'] sage: K.<T> = Frac(A) sage: phi = DrinfeldModule(A, [T, 1]) sage: phi.goss_polynomial(1) X sage: phi.goss_polynomial(2) X^2 sage: phi.goss_polynomial(3) X^3 + (1/(T^2 + T))*X^2 - >>> from sage.all import * >>> A = GF(Integer(2))['T'] >>> K = Frac(A, names=('T',)); (T,) = K._first_ngens(1) >>> phi = DrinfeldModule(A, [T, Integer(1)]) >>> phi.goss_polynomial(Integer(1)) X >>> phi.goss_polynomial(Integer(2)) X^2 >>> phi.goss_polynomial(Integer(3)) X^3 + (1/(T^2 + T))*X^2 - Base fields of \(\mathbb{F}_q[T]\)-characteristic zero - The base fields need not only be fraction fields of polynomials ring. In the following example, we construct a Drinfeld module over \(\mathbb{F}_q((1/T))\), the completion of the rational function field at the place \(1/T\): - sage: A.<T> = GF(2)[] sage: L.<s> = LaurentSeriesRing(GF(2)) # s = 1/T sage: phi = DrinfeldModule(A, [1/s, s + s^2 + s^5 + O(s^6), 1+1/s]) sage: phi(T) (s^-1 + 1)*t^2 + (s + s^2 + s^5 + O(s^6))*t + s^-1 - >>> from sage.all import * >>> A = GF(Integer(2))['T']; (T,) = A._first_ngens(1) >>> L = LaurentSeriesRing(GF(Integer(2)), names=('s',)); (s,) = L._first_ngens(1)# s = 1/T >>> phi = DrinfeldModule(A, [Integer(1)/s, s + s**Integer(2) + s**Integer(5) + O(s**Integer(6)), Integer(1)+Integer(1)/s]) >>> phi(T) (s^-1 + 1)*t^2 + (s + s^2 + s^5 + O(s^6))*t + s^-1 - One can also construct Drinfeld modules over SageMath’s global function fields: - sage: A.<T> = GF(5)[] sage: K.<z> = FunctionField(GF(5)) # z = T sage: phi = DrinfeldModule(A, [z, 1, z^2]) sage: phi(T) z^2*t^2 + t + z - >>> from sage.all import * >>> A = GF(Integer(5))['T']; (T,) = A._first_ngens(1) >>> K = FunctionField(GF(Integer(5)), names=('z',)); (z,) = K._first_ngens(1)# z = T >>> phi = DrinfeldModule(A, [z, Integer(1), z**Integer(2)]) >>> phi(T) z^2*t^2 + t + z - exponential(prec=+Infinity, name='z')[source]¶
- Return the exponential of this Drinfeld module. - Note that the exponential is only defined when the \(\mathbb{F}_q[T]\)-characteristic is zero. - INPUT: - prec– an integer or- Infinity(default:- Infinity); the precision at which the series is returned; if- Infinity, a lazy power series in returned, else, a classical power series is returned.
- name– string (default:- 'z'); the name of the generator of the lazy power series ring
 - EXAMPLES: - sage: A = GF(2)['T'] sage: K.<T> = Frac(A) sage: phi = DrinfeldModule(A, [T, 1]) sage: q = A.base_ring().cardinality() - >>> from sage.all import * >>> A = GF(Integer(2))['T'] >>> K = Frac(A, names=('T',)); (T,) = K._first_ngens(1) >>> phi = DrinfeldModule(A, [T, Integer(1)]) >>> q = A.base_ring().cardinality() - When - precis- Infinity(which is the default), the exponential is returned as a lazy power series, meaning that any of its coefficients can be computed on demands:- sage: exp = phi.exponential(); exp z + ((1/(T^2+T))*z^2) + ((1/(T^8+T^6+T^5+T^3))*z^4) + O(z^8) sage: exp[2^4] 1/(T^64 + T^56 + T^52 + ... + T^27 + T^23 + T^15) sage: exp[2^5] 1/(T^160 + T^144 + T^136 + ... + T^55 + T^47 + T^31) - >>> from sage.all import * >>> exp = phi.exponential(); exp z + ((1/(T^2+T))*z^2) + ((1/(T^8+T^6+T^5+T^3))*z^4) + O(z^8) >>> exp[Integer(2)**Integer(4)] 1/(T^64 + T^56 + T^52 + ... + T^27 + T^23 + T^15) >>> exp[Integer(2)**Integer(5)] 1/(T^160 + T^144 + T^136 + ... + T^55 + T^47 + T^31) - On the contrary, when - precis a finite number, all the required coefficients are computed at once:- sage: phi.exponential(prec=10) z + (1/(T^2 + T))*z^2 + (1/(T^8 + T^6 + T^5 + T^3))*z^4 + (1/(T^24 + T^20 + T^18 + T^17 + T^14 + T^13 + T^11 + T^7))*z^8 + O(z^10) - >>> from sage.all import * >>> phi.exponential(prec=Integer(10)) z + (1/(T^2 + T))*z^2 + (1/(T^8 + T^6 + T^5 + T^3))*z^4 + (1/(T^24 + T^20 + T^18 + T^17 + T^14 + T^13 + T^11 + T^7))*z^8 + O(z^10) - Example in higher rank: - sage: A = GF(5)['T'] sage: K.<T> = Frac(A) sage: phi = DrinfeldModule(A, [T, T^2, T + T^2 + T^4, 1]) sage: exp = phi.exponential(); exp z + ((T/(T^4+4))*z^5) + O(z^8) - >>> from sage.all import * >>> A = GF(Integer(5))['T'] >>> K = Frac(A, names=('T',)); (T,) = K._first_ngens(1) >>> phi = DrinfeldModule(A, [T, T**Integer(2), T + T**Integer(2) + T**Integer(4), Integer(1)]) >>> exp = phi.exponential(); exp z + ((T/(T^4+4))*z^5) + O(z^8) - The exponential is the compositional inverse of the logarithm (see - logarithm()):- sage: log = phi.logarithm(); log z + ((4*T/(T^4+4))*z^5) + O(z^8) sage: exp.compose(log) z + O(z^8) sage: log.compose(exp) z + O(z^8) - >>> from sage.all import * >>> log = phi.logarithm(); log z + ((4*T/(T^4+4))*z^5) + O(z^8) >>> exp.compose(log) z + O(z^8) >>> log.compose(exp) z + O(z^8) - REFERENCE: - See section 4.6 of [Gos1998] for the definition of the exponential. 
 - goss_polynomial(n, var='X')[source]¶
- 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– string (default:- 'X'); the name of polynomial variable
 - OUTPUT: a univariate polynomial in - varover 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 - >>> from sage.all import * >>> A = GF(Integer(3))['T'] >>> K = Frac(A, names=('T',)); (T,) = K._first_ngens(1) >>> phi = DrinfeldModule(A, [T, Integer(1)]) # Carlitz module >>> phi.goss_polynomial(Integer(1)) X >>> phi.goss_polynomial(Integer(2)) X^2 >>> phi.goss_polynomial(Integer(4)) X^4 + (1/(T^3 + 2*T))*X^2 >>> phi.goss_polynomial(Integer(5)) X^5 + (2/(T^3 + 2*T))*X^3 >>> phi.goss_polynomial(Integer(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 - REFERENCE: - Section 3 of [Gek1988] provides an exposition of Goss polynomials. 
 - logarithm(prec=+Infinity, name='z')[source]¶
- Return the logarithm of the given Drinfeld module. - By definition, the logarithm is the compositional inverse of the exponential (see - exponential()). Note that the logarithm is only defined when the \(\mathbb{F}_q[T]\)-characteristic is zero.- INPUT: - prec– an integer or- Infinity(default:- Infinity); the precision at which the series is returned; if- Infinity, a lazy power series in returned
- name– string (default:- 'z'); the name of the generator of the lazy power series ring
 - EXAMPLES: - sage: A = GF(2)['T'] sage: K.<T> = Frac(A) sage: phi = DrinfeldModule(A, [T, 1]) - >>> from sage.all import * >>> A = GF(Integer(2))['T'] >>> K = Frac(A, names=('T',)); (T,) = K._first_ngens(1) >>> phi = DrinfeldModule(A, [T, Integer(1)]) - When - precis- Infinity(which is the default), the logarithm is returned as a lazy power series, meaning that any of its coefficients can be computed on demands:- sage: log = phi.logarithm(); log z + ((1/(T^2+T))*z^2) + ((1/(T^6+T^5+T^3+T^2))*z^4) + O(z^8) sage: log[2^4] 1/(T^30 + T^29 + T^27 + ... + T^7 + T^5 + T^4) sage: log[2^5] 1/(T^62 + T^61 + T^59 + ... + T^8 + T^6 + T^5) - >>> from sage.all import * >>> log = phi.logarithm(); log z + ((1/(T^2+T))*z^2) + ((1/(T^6+T^5+T^3+T^2))*z^4) + O(z^8) >>> log[Integer(2)**Integer(4)] 1/(T^30 + T^29 + T^27 + ... + T^7 + T^5 + T^4) >>> log[Integer(2)**Integer(5)] 1/(T^62 + T^61 + T^59 + ... + T^8 + T^6 + T^5) - If - precis a finite number, all the required coefficients are computed at once:- sage: phi.logarithm(prec=10) z + (1/(T^2 + T))*z^2 + (1/(T^6 + T^5 + T^3 + T^2))*z^4 + (1/(T^14 + T^13 + T^11 + T^10 + T^7 + T^6 + T^4 + T^3))*z^8 + O(z^10) - >>> from sage.all import * >>> phi.logarithm(prec=Integer(10)) z + (1/(T^2 + T))*z^2 + (1/(T^6 + T^5 + T^3 + T^2))*z^4 + (1/(T^14 + T^13 + T^11 + T^10 + T^7 + T^6 + T^4 + T^3))*z^8 + O(z^10) - Example in higher rank: - sage: A = GF(5)['T'] sage: K.<T> = Frac(A) sage: phi = DrinfeldModule(A, [T, T^2, T + T^2 + T^4, 1]) sage: phi.logarithm() z + ((4*T/(T^4+4))*z^5) + O(z^8) - >>> from sage.all import * >>> A = GF(Integer(5))['T'] >>> K = Frac(A, names=('T',)); (T,) = K._first_ngens(1) >>> phi = DrinfeldModule(A, [T, T**Integer(2), T + T**Integer(2) + T**Integer(4), Integer(1)]) >>> phi.logarithm() z + ((4*T/(T^4+4))*z^5) + O(z^8) 
 
- class sage.rings.function_field.drinfeld_modules.charzero_drinfeld_module.DrinfeldModule_rational(gen, category)[source]¶
- Bases: - DrinfeldModule_charzero- A class for Drinfeld modules defined over the fraction field of the underlying function field. - class_polynomial()[source]¶
- Return the class polynomial, that is the Fitting ideal of the class module, of this Drinfeld module. - We refer to [Tae2012] for the definition and basic properties of the class module. - EXAMPLES: - We check that the class module of the Carlitz module is trivial: - sage: q = 5 sage: Fq = GF(q) sage: A = Fq['T'] sage: K.<T> = Frac(A) sage: C = DrinfeldModule(A, [T, 1]); C Drinfeld module defined by T |--> t + T sage: C.class_polynomial() 1 - >>> from sage.all import * >>> q = Integer(5) >>> Fq = GF(q) >>> A = Fq['T'] >>> K = Frac(A, names=('T',)); (T,) = K._first_ngens(1) >>> C = DrinfeldModule(A, [T, Integer(1)]); C Drinfeld module defined by T |--> t + T >>> C.class_polynomial() 1 - When the coefficients of the Drinfeld module have small enough degrees, the class module is always trivial: - sage: gs = [T] + [A.random_element(degree = q^i) ....: for i in range(1, 5)] sage: phi = DrinfeldModule(A, gs) sage: phi.class_polynomial() 1 - >>> from sage.all import * >>> gs = [T] + [A.random_element(degree = q**i) ... for i in range(Integer(1), Integer(5))] >>> phi = DrinfeldModule(A, gs) >>> phi.class_polynomial() 1 - Here is an example with a nontrivial class module: - sage: phi = DrinfeldModule(A, [T, 2*T^14 + 2*T^4]) sage: phi.class_polynomial() T + 3 - >>> from sage.all import * >>> phi = DrinfeldModule(A, [T, Integer(2)*T**Integer(14) + Integer(2)*T**Integer(4)]) >>> phi.class_polynomial() T + 3 
 - coefficient_in_function_ring(n)[source]¶
- Return the \(n\)-th coefficient of this Drinfeld module as an element of the underlying function ring. - INPUT: - n– an integer
 - EXAMPLES: - sage: q = 5 sage: Fq = GF(q) sage: A = Fq['T'] sage: R = Fq['U'] sage: K.<U> = Frac(R) sage: phi = DrinfeldModule(A, [U, 0, U^2, U^3]) sage: phi.coefficient_in_function_ring(2) T^2 - >>> from sage.all import * >>> q = Integer(5) >>> Fq = GF(q) >>> A = Fq['T'] >>> R = Fq['U'] >>> K = Frac(R, names=('U',)); (U,) = K._first_ngens(1) >>> phi = DrinfeldModule(A, [U, Integer(0), U**Integer(2), U**Integer(3)]) >>> phi.coefficient_in_function_ring(Integer(2)) T^2 - Compare with the method meth:\(coefficient\): - sage: phi.coefficient(2) U^2 - >>> from sage.all import * >>> phi.coefficient(Integer(2)) U^2 - If the required coefficient is not a polynomials, an error is raised: - sage: psi = DrinfeldModule(A, [U, 1/U]) sage: psi.coefficient_in_function_ring(0) T sage: psi.coefficient_in_function_ring(1) Traceback (most recent call last): ... ValueError: coefficient is not polynomial - >>> from sage.all import * >>> psi = DrinfeldModule(A, [U, Integer(1)/U]) >>> psi.coefficient_in_function_ring(Integer(0)) T >>> psi.coefficient_in_function_ring(Integer(1)) Traceback (most recent call last): ... ValueError: coefficient is not polynomial 
 - coefficients_in_function_ring(sparse=True)[source]¶
- Return the coefficients of this Drinfeld module as elements of the underlying function ring. - INPUT: - sparse– a boolean (default:- True); if- True, only return the nonzero coefficients; otherwise, return all of them.
 - EXAMPLES: - sage: q = 5 sage: Fq = GF(q) sage: A = Fq['T'] sage: R = Fq['U'] sage: K.<U> = Frac(R) sage: phi = DrinfeldModule(A, [U, 0, U^2, U^3]) sage: phi.coefficients_in_function_ring() [T, T^2, T^3] sage: phi.coefficients_in_function_ring(sparse=False) [T, 0, T^2, T^3] - >>> from sage.all import * >>> q = Integer(5) >>> Fq = GF(q) >>> A = Fq['T'] >>> R = Fq['U'] >>> K = Frac(R, names=('U',)); (U,) = K._first_ngens(1) >>> phi = DrinfeldModule(A, [U, Integer(0), U**Integer(2), U**Integer(3)]) >>> phi.coefficients_in_function_ring() [T, T^2, T^3] >>> phi.coefficients_in_function_ring(sparse=False) [T, 0, T^2, T^3] - Compare with the method meth:\(coefficients\): - sage: phi.coefficients() [U, U^2, U^3] - >>> from sage.all import * >>> phi.coefficients() [U, U^2, U^3] - If the coefficients are not polynomials, an error is raised: - sage: psi = DrinfeldModule(A, [U, 1/U]) sage: psi.coefficients_in_function_ring() Traceback (most recent call last): ... ValueError: coefficients are not polynomials - >>> from sage.all import * >>> psi = DrinfeldModule(A, [U, Integer(1)/U]) >>> psi.coefficients_in_function_ring() Traceback (most recent call last): ... ValueError: coefficients are not polynomials