Implemented pseudo_hermite_form

This commit is contained in:
Montessinos Mickael Gerard Bernard 2024-03-01 17:07:42 +02:00
parent 3c37da73b0
commit 783f1d6371
1 changed files with 190 additions and 2 deletions

View File

@ -1,3 +1,13 @@
r"""
Implementations of all sort of algorithms for function field that are not
yet implemented in Sage.
REFERENCE:
.. [Coh00] H. Cohen
Advanced topics in computational number theory
Springer
2000
"""
###########################################################################
# Copyright (C) 2024 Mickaël Montessinos (mickael.montessinos@mif.vu.lt),#
# #
@ -19,6 +29,7 @@ from sage.rings.function_field.function_field_rational\
import RationalFunctionField
from sage.rings.function_field.order_rational\
import FunctionFieldMaximalOrderInfinite_rational
from sage.rings.function_field.order import FunctionFieldOrderInfinite
@cached_function
def all_infinite_places(K):
@ -116,11 +127,16 @@ def infinite_hermite_form(mat,include_zero_cols=True,transformation=False):
sage: K.<x> = FunctionField(GF(3))
sage: R = K.maximal_order_infinite()
sage: mat = matrix(R,[[1, x**-1, x**-2, (x**3+1) / x**3], [(2*x+2) / (x**3+2), x**-2, (x**2+2) / (x**4+1), 1]])
sage: H,T = infinite_hermite_form(mat,transformation=True); H
sage: H, T = infinite_hermite_form(mat,transformation=True); H
[0 0 1 0]
[0 0 0 1]
sage: mat*T == H
True
sage: H, T = infinite_hermite_form(mat, False, True); H
[1 0]
[0 1]
sage: mat*T == H
True
TESTS:
@ -165,6 +181,7 @@ def infinite_hermite_form(mat,include_zero_cols=True,transformation=False):
H *= E
if not include_zero_cols:
H = H[:,r-n:]
T = T[:,r-n:]
if transformation:
return H,T
return H
@ -241,7 +258,7 @@ def infinite_order_xgcd(ideals):
H,U = infinite_hermite_form(C, include_zero_cols=False, transformation=True)
if not (H/den).is_one():
raise ValueError("The ideals should be coprime.")
v = U[:,-n].list()
v = U[:,0].list()
return [sum([ideals_bases[i][j]*v[n*i+j] for j in range(n)]) for i in range(k)]
@ -464,3 +481,174 @@ def smallest_norm_first(mat,i = 0,norms=[]):
norms[i] = norms[j+i]
norms[j+i] = n
return norms
def finite_order_xgcd(left, right):
r"""
Compute `a \in left` and `b \in right` such that `a + b = 1`
INPUT:
-a: FunctionFieldIdeal_polymod: ideal of a finite maximal order
-b: FunctionFieldIdeal_polymod: ideal of the same maximal order
ALGORITHM:
[Coh00]_ Algorithm 1.3.2
EXAMPLES ::
sage: from vector_bundle.function_field_utility import finite_order_xgcd
sage: F.<x> = FunctionField(GF(7))
sage: R.<y> = F[]
sage: K.<y> = F.extension(y^2 - x^3 - x)
sage: places = K.places_finite()
sage: left = places[0].prime_ideal()
sage: right = places[1].prime_ideal()
sage: a, b = finite_order_xgcd(left, right)
sage: a in left
True
sage: b in right
True
sage: a + b
1
"""
O = left.base_ring()
O_basis = O.basis()
m_left = left.hnf()
n = m_left.nrows()
m_right = right.hnf()
c = block_matrix([[m_left], [m_right]])
h, u = c.hermite_form(False, True)
if not h.is_one():
raise ValueError('The ideals left and right must be coprime.')
x = u[0,:n].list()
a = sum([c*e for c, e in zip(x,left.gens_over_base())])
return a, 1-a
def euclidean_step(ideal_a, ideal_b, a, b, d=None):
r"""
Let `d = ideal_a a + ideal_b b`. Return `u \in ideal_a d^-1` and
`v \in ideal_b d^-1` such that au + bv = 1
ALGORITHM:
[Coh00]_ Theorem 1.3.3
EXAMPLES ::
sage: from vector_bundle.function_field_utility import euclidean_step
sage: F.<x> = FunctionField(GF(7))
sage: R.<y> = F[]
sage: K.<y> = F.extension(y^2 - x^3 - x)
sage: ideals = [P.prime_ideal() for P in K.places_finite()[:2]]
sage: a = x^2*y + 3
sage: b = 2*y*x^5
sage: d = a*ideals[0] + b*ideals[1]
sage: u, v = euclidean_step(ideals[0], ideals[1], a, b)
sage: u in ideals[0] * d^-1
True
sage: v in ideals[1] * d^-1
True
sage: a*u + b*v
1
"""
infinite = isinstance(ideal_a.base_ring(),FunctionFieldOrderInfinite)
if a == 0:
return 0, b^-1
if b == 0:
return a^-1, 0
if d == None:
d = a*ideal_a + b*ideal_b
I = a * ideal_a * d**-1
J = b * ideal_b * d**-1
#It would make sense to distinguish between finite and infinite order
#in a unified xgcd function, but the hnf form for infinite ideals
#needs to be refactored: we currently use opposite convention from
#that of sage.
if infinite:
s, t = infinite_order_xgcd([I,J])
else:
s, t = finite_order_xgcd(I, J)
return s/a, t/b
def pseudo_hermite_form(ideals, mat, include_zero_cols=True, transformation=False):
r"""
Return the hermite form of the pseudo-matrix ``(ideals, mat)`` with
coefficients in a function field and ideals in a maximal order.
If the maximal order is infinite, the ideals of the pseudo-matrix output
are all trivial.
WARNING:
Uses the opposite convention from sage for hermite forms, aligns with
Cohen's book instead.
The final reduction step is not implemented as having triangular matrices
is enough for our purposes.
ALGORITHM:
- Algorithm 1.4.7 from [Coh00]_
EXAMPLES ::
sage: from vector_bundle.function_field_utility import pseudo_hermite_form
sage: F.<x> = FunctionField(GF(7))
sage: R.<y> = F[]
sage: K.<y> = F.extension(y^2 - x^3 - x)
sage: ideals = [P.prime_ideal() for P in K.places_finite()[:3]]
sage: mat = matrix(K, [[1, x, y],[2, x+1, y+1]])
sage: h_ideals, h, u = pseudo_hermite_form(ideals, mat, transformation=True)
sage: h
[ 0 1 (2*x + 5)*y + 4*x]
[ 0 0 1]
sage: h == mat * u
True
sage: all([u[i,j] in ideals[i] * h_ideals[j]^-1 for i in range(3) for j in range(3)])
True
sage: prod(ideals) == u.determinant() * prod(h_ideals)
True
"""
K = mat.base_ring()
k = mat.ncols()
n = mat.nrows()
U = identity_matrix(K,k)
h = copy(mat)
h_ideals = copy(ideals)
j = k-1
for i in range(n-1, -1, -1):
#Check zero
m = [h[i,m] == 0 for m in range(j, -1, -1)].index(False)
h[:,m], h[:,j] = h[:,j], h[:,m]
U[:,m], U[:,j] = U[:,j], U[:,m]
h_ideals[m], h_ideals[j] = h_ideals[j], h_ideals[m]
#Put 1 on the main diagonal
a = h[i,j]**-1
h_ideals[j] *= h[i,j]
h[:,j] *= a
U[:,j] *= a
for m in range(j-1,-1,-1):
if h[i,m] == 0:
continue
#Euclidean step
partial = h[i,m]*h_ideals[m] + h_ideals[j]
u, v = euclidean_step(h_ideals[m], h_ideals[j],
h[i,m], 1, partial)
U[:, m], U[:, j] = U[:, m] - h[i, m]*U[:, j], u*U[:, m] + v*U[:, j]
h[:, m], h[:, j] = h[:, m] - h[i, m]*h[:, j], u*h[:, m] + v*h[:, j]
h_ideals[m], h_ideals[j] = h_ideals[m] * h_ideals[j] * partial**-1, partial
#Final reduction of row i
#Not implemented for the moment since we only need triangular matrices
#for our current purposes.
j -=1
if not include_zero_cols:
h = h[:,k-n:]
h_ideals = h_ideals[k-n:]
U = U[:,k-n:]
if transformation:
return h_ideals, h, U
return h_ideals, h