414 lines
14 KiB
Python
414 lines
14 KiB
Python
class superelliptic_witt:
|
|
def __init__(self, t, f):
|
|
''' Define Witt function on C of the form [t] + V(f). '''
|
|
self.curve = t.curve
|
|
C = t.curve
|
|
p = C.characteristic
|
|
self.t = t #superelliptic_function
|
|
self.f = f #superelliptic_function
|
|
|
|
def __repr__(self):
|
|
f = self.f
|
|
t = self.t
|
|
if f.function == 0:
|
|
return "[" + str(t) + "]"
|
|
if t.function == 0:
|
|
return "V(" + str(f) + ")"
|
|
return "[" + str(t) + "] + V(" + str(f) + ")"
|
|
|
|
def __neg__(self):
|
|
f = self.f
|
|
t = self.t
|
|
return superelliptic_witt(-t, -f)
|
|
|
|
def __add__(self, other):
|
|
C = self.curve
|
|
second_coor = 0*C.x
|
|
X = self.t
|
|
Y = other.t
|
|
for i in range(1, p):
|
|
second_coor -= binomial_prim(p, i)*X^i*Y^(p-i)
|
|
return superelliptic_witt(self.t + other.t, self.f + other.f + second_coor)
|
|
|
|
def __sub__(self, other):
|
|
return self + (-other)
|
|
|
|
def __rmul__(self, other):
|
|
p = self.curve.characteristic
|
|
if other in ZZ:
|
|
if other == 0:
|
|
return superelliptic_witt(0*C.x, 0*C.x)
|
|
if other > 0:
|
|
return self + (other-1)*self
|
|
if other < 0:
|
|
return (-other)*(-self)
|
|
if other in QQ:
|
|
other_integer = Integers(p^2)(other)
|
|
return other_integer*self
|
|
|
|
def __mul__(self, other):
|
|
C = self.curve
|
|
p = C.characteristic
|
|
if isinstance(other, superelliptic_witt):
|
|
t1 = self.t
|
|
f1 = self.f
|
|
t2 = other.t
|
|
f2 = other.f
|
|
return superelliptic_witt(t1*t2, t1^p*f2 + t2^p*f1)
|
|
if isinstance(other, superelliptic_drw_form):
|
|
h1 = other.h1
|
|
h2 = other.h2
|
|
omega = other.omega
|
|
t = self.t
|
|
f = self.f
|
|
aux_form = t^p*omega - h2*t^(p-1)*t.diffn() + f*h1^p*(C.x)^(p-1)*C.dx
|
|
return superelliptic_drw_form(t*h1, aux_form, t^p*h2)
|
|
|
|
def __eq__(self, other):
|
|
return self.t == other.t and self.f == other.f
|
|
|
|
def diffn(self, dy_w = 0):
|
|
if dy_w == 0:
|
|
dy_w = self.curve.dy_w()
|
|
C = self.curve
|
|
t = self.t
|
|
f = self.f
|
|
fC = C.polynomial
|
|
F = C.base_ring
|
|
Rxy.<x, y> = PolynomialRing(F, 2)
|
|
if t.function == 0:
|
|
return superelliptic_drw_form(0*C.x, 0*C.dx, f)
|
|
t_polynomial = t.function
|
|
num = t_polynomial.numerator()
|
|
den = t_polynomial.denominator()
|
|
num_t_fct = superelliptic_function(C, num)
|
|
den_t_fct = superelliptic_function(C, den)
|
|
inv_den_t_fct = superelliptic_function(C, 1/den)
|
|
if den != 1:
|
|
# d([N/D] + V(f)) = [1/D]*d([N]) - [N]*[D^(-2)]*d([D]) + dV(f)
|
|
return ((den_t_fct)^(-1)).teichmuller()*num_t_fct.teichmuller().diffn() - ((den_t_fct)^(-2)).teichmuller()*num_t_fct.teichmuller()*den_t_fct.teichmuller().diffn() + superelliptic_drw_form(0*C.x, 0*C.dx, f)
|
|
t_polynomial = Rxy(t_polynomial)
|
|
M = t_polynomial.monomials()[0]
|
|
a = t_polynomial.monomial_coefficient(M)
|
|
#[P] = [aM] + Q, where Q = ([P] - [aM])
|
|
aM_fct = superelliptic_function(C, a*M)
|
|
Q = self - aM_fct.teichmuller()
|
|
exp_x = M.exponents()[0][0]
|
|
exp_y = M.exponents()[0][1]
|
|
return Q.diffn() + exp_x*superelliptic_drw_form(aM_fct/C.x, 0*C.dx, 0*C.x) + exp_y*(aM_fct/C.y).teichmuller()*dy_w
|
|
|
|
def binomial_prim(p, i):
|
|
return binomial(p, i)/p
|
|
|
|
def reduce_rational_fct(fct, p):
|
|
Rxy.<x, y> = PolynomialRing(QQ)
|
|
Fxy = FractionField(Rxy)
|
|
fct = Fxy(fct)
|
|
num = Rxy(fct.numerator())
|
|
den = Rxy(fct.denominator())
|
|
num1 = Rxy(0)
|
|
for m in num.monomials():
|
|
a = num.monomial_coefficient(m)
|
|
num1 += (a%p^2)*m
|
|
den1 = Rxy(0)
|
|
for m in den.monomials():
|
|
a = den.monomial_coefficient(m)
|
|
den1 += (a%p^2)*m
|
|
return num1/den1
|
|
|
|
def teichmuller(fct):
|
|
C = fct.curve
|
|
return superelliptic_witt(fct, 0*C.x)
|
|
|
|
superelliptic_function.teichmuller = teichmuller
|
|
|
|
#dy = [f(x)]'/2*y dx
|
|
#[f1 + M] = [f1] + [M] + V(cos)
|
|
#d[f1 + M] = d[f1] + d[M] + dV(f1*M)
|
|
#M = b x^a
|
|
#d[M] = a*[b x^(a-1)]
|
|
|
|
def auxilliary_derivative(P):
|
|
'''Return "derivative" of P, where P depends only on x. In other words d[P(x)].'''
|
|
P0 = P.t.function
|
|
P1 = P.f.function
|
|
C = P.curve
|
|
F = C.base_ring
|
|
Rx.<x> = PolynomialRing(F)
|
|
P0 = Rx(P0)
|
|
P1 = Rx(P1)
|
|
if P0 == 0:
|
|
return superelliptic_drw_form(0*C.x, 0*C.dx, P.f)
|
|
M = P0.monomials()[0]
|
|
a = P0.monomial_coefficient(M)
|
|
#[P] = [aM] + Q, where Q = ([P] - [aM])
|
|
aM_fct = superelliptic_function(C, a*M)
|
|
Q = P - aM_fct.teichmuller()
|
|
exp = M.exponents()[0]
|
|
return auxilliary_derivative(Q) + exp*superelliptic_drw_form(aM_fct/C.x, 0*C.dx, 0*C.x)
|
|
|
|
class superelliptic_drw_form:
|
|
def __init__(self, h1, omega, h2):
|
|
'''Form [h1] d[x] + V(omega) + dV([h])'''
|
|
self.curve = h1.curve
|
|
self.h1 = h1
|
|
self.omega = omega
|
|
self.h2 = h2
|
|
|
|
def r(self):
|
|
C = self.curve
|
|
h1 = self.h1
|
|
return superelliptic_form(C, h1.function)
|
|
|
|
def __eq__(self, other):
|
|
eq1 = (self.h1 == self.h1)
|
|
try:
|
|
H = (self.h2 - other.h2).pth_root()
|
|
except:
|
|
return False
|
|
eq2 = ((self.omega - other.omega).cartier() - H.diffn()) == 0*self.curve.dx
|
|
if eq1 and eq2:
|
|
return True
|
|
return False
|
|
|
|
def __repr__(self):
|
|
C = self.curve
|
|
h1 = self.h1
|
|
omega = self.omega
|
|
h2 = self.h2
|
|
result = ""
|
|
if h1.function != 0:
|
|
result += "[" + str(h1) + "] d[x]"
|
|
if h1.function !=0 and omega.form != 0:
|
|
result += " + "
|
|
if omega.form != 0:
|
|
result += "V(" + str(omega) + ")"
|
|
if h2.function !=0 and omega.form != 0:
|
|
result += " + "
|
|
if h2.function != 0:
|
|
result += "dV([" + str(h2) +"])"
|
|
if result == "":
|
|
result += "0"
|
|
return result
|
|
|
|
def __rmul__(self, other):
|
|
h1 = self.h1
|
|
h2 = self.h2
|
|
omega = self.omega
|
|
p = self.curve.characteristic
|
|
if other in ZZ:
|
|
if other == 0:
|
|
return superelliptic_drw_form(0*C.x, 0*C.dx, 0*C.x)
|
|
if other > 0:
|
|
return self + (other-1)*self
|
|
if other < 0:
|
|
return (-other)*(-self)
|
|
if other in QQ:
|
|
other_integer = Integers(p^2)(other)
|
|
return other_integer*self
|
|
t = other.t
|
|
f = other.f
|
|
aux_form = t^p*omega - h2*t^(p-1)*t.diffn() + f*h1^p*(C.x)^(p-1)*C.dx
|
|
return superelliptic_drw_form(t*h1, aux_form, t^p*h2)
|
|
|
|
def __neg__(self):
|
|
C = self.curve
|
|
h1 = self.h1
|
|
h2 = self.h2
|
|
omega = self.omega
|
|
return superelliptic_drw_form(-h1, -omega, -h2)
|
|
|
|
def __sub__(self, other):
|
|
return self + (-other)
|
|
|
|
def __add__(self, other):
|
|
C = self.curve
|
|
h1 = self.h1
|
|
h2 = self.h2
|
|
omega = self.omega
|
|
H1 = other.h1
|
|
H2 = other.h2
|
|
OMEGA = other.omega
|
|
aux = (teichmuller(h1) + teichmuller(H1))*superelliptic_drw_form(C.one, 0*C.dx, 0*C.x)
|
|
h1_aux = aux.h1
|
|
h2_aux = aux.h2
|
|
omega_aux = aux.omega
|
|
return superelliptic_drw_form(h1_aux, omega + OMEGA + omega_aux, h2 + H2 + h2_aux)
|
|
|
|
def frobenius(self):
|
|
C = self.curve
|
|
h1 = self.h1
|
|
h2 = self.h2
|
|
p = C.characteristic
|
|
return h1^p*C.x^(p-1)*C.dx + h2.diffn()
|
|
|
|
def mult_by_p(omega):
|
|
C = omega.curve
|
|
fct = omega.form
|
|
Fxy, Rxy, x, y = C.fct_field
|
|
omega = superelliptic_form(C, fct^p * x^(p-1))
|
|
result = superelliptic_drw_form(0*C.x, omega, 0*C.x)
|
|
return result
|
|
|
|
def verschiebung(elt):
|
|
C = elt.curve
|
|
if isinstance(elt, superelliptic_function):
|
|
return superelliptic_witt(0*C.x, elt)
|
|
|
|
if isinstance(elt, superelliptic_form):
|
|
return superelliptic_drw_form(0*C.x, elt, 0*C.x)
|
|
|
|
superelliptic_form.verschiebung = verschiebung
|
|
superelliptic_function.verschiebung = verschiebung
|
|
|
|
|
|
class superelliptic_drw_cech:
|
|
def __init__(self, omega0, f):
|
|
self.curve = omega0.curve
|
|
self.omega0 = omega0
|
|
self.omega8 = omega0 - f.diffn()
|
|
self.f = f
|
|
|
|
|
|
def reduce(self):
|
|
C = self.curve
|
|
fct = self.f
|
|
f_first_comp = fct.t
|
|
f_second_comp = fct.f
|
|
decomp_first_comp = decomposition_g0_g8(f_first_comp)
|
|
decomp_second_comp = decomposition_g0_g8(f_second_comp)
|
|
new = self
|
|
new.omega0 -= decomposition_g0_g8(f_first_comp)[0].teichmuller().diffn()
|
|
new.omega0 -= decomposition_g0_g8(f_second_comp)[0].verschiebung().diffn()
|
|
new.f = decomposition_g0_g8(f_first_comp)[2].teichmuller() + decomposition_g0_g8(f_second_comp)[2].verschiebung()
|
|
new.omega8 = new.omega0 - new.f.diffn()
|
|
return new
|
|
|
|
def __repr__(self):
|
|
return("(" + str(self.omega0) + ", "+ str(self.f) + ", " + str(self.omega8) + ")")
|
|
|
|
def __add__(self, other):
|
|
C = self.curve
|
|
omega0 = self.omega0
|
|
f = self.f
|
|
omega0_1 = other.omega0
|
|
f_1 = other.f
|
|
return superelliptic_drw_cech(omega0 + omega0_1, f + f_1)
|
|
|
|
def __sub__(self, other):
|
|
C = self.curve
|
|
omega0 = self.omega0
|
|
f = self.f
|
|
omega0_1 = other.omega0
|
|
f_1 = other.f
|
|
return superelliptic_drw_cech(omega0 - omega0_1, f - f_1)
|
|
|
|
def __neg__(self):
|
|
C = self.curve
|
|
omega0 = self.omega0
|
|
f = self.f
|
|
return superelliptic_drw_cech(-omega0, -f)
|
|
|
|
def __rmul__(self, other):
|
|
omega0 = self.omega0
|
|
f = self.f
|
|
return superelliptic_drw_cech(other*omega0, other*f)
|
|
|
|
def r(self):
|
|
omega0 = self.omega0
|
|
f = self.f
|
|
C = self.curve
|
|
return superelliptic_cech(C, omega0.h1*C.dx, f.t)
|
|
|
|
def coordinates(self, basis = 0):
|
|
C = self.curve
|
|
g = C.genus()
|
|
coord_mod_p = self.r().coordinates()
|
|
print(coord_mod_p)
|
|
coord_lifted = [lift(a) for a in coord_mod_p]
|
|
if basis == 0:
|
|
basis = C.crystalline_cohomology_basis()
|
|
aux = self
|
|
for i, a in enumerate(basis):
|
|
aux -= coord_lifted[i]*a
|
|
print('aux before reduce', aux)
|
|
#aux = aux.reduce() # Now aux = p*cech class.
|
|
# Now aux should be of the form (V(smth), V(smth), V(smth))
|
|
print('aux V(smth)', aux)
|
|
aux_divided_by_p = superelliptic_cech(C, aux.omega0.omega.cartier(), aux.f.f.pth_root())
|
|
print('aux.omega0.omega.cartier()', aux.omega0.omega.cartier())
|
|
coord_aux_divided_by_p = aux_divided_by_p.coordinates()
|
|
coord_aux_divided_by_p = [ZZ(a) for a in coord_aux_divided_by_p]
|
|
coordinates = [ (coord_lifted[i] + p*coord_aux_divided_by_p[i])%p^2 for i in range(2*g)]
|
|
return coordinates
|
|
|
|
def is_regular(self):
|
|
print(self.omega0.r().is_regular_on_U0(), self.omega8.r().is_regular_on_Uinfty(), self.omega0.frobenius().is_regular_on_U0(), self.omega8.frobenius().is_regular_on_Uinfty())
|
|
eq1 = self.omega0.r().is_regular_on_U0() and self.omega8.r().is_regular_on_Uinfty()
|
|
eq2 = self.omega0.frobenius().is_regular_on_U0() and self.omega8.frobenius().is_regular_on_Uinfty()
|
|
return eq1 and eq2
|
|
|
|
|
|
def de_rham_witt_lift(cech_class, prec = 50):
|
|
C = cech_class.curve
|
|
g = C.genus()
|
|
omega0 = cech_class.omega0
|
|
omega8 = cech_class.omega8
|
|
fct = cech_class.f
|
|
omega0_regular = regular_form(omega0) #Present omega0 in the form P dx + Q dy
|
|
omega0_lift = omega0_regular[0].teichmuller()*(C.x.teichmuller().diffn()) + omega0_regular[1].teichmuller()*(C.y.teichmuller().diffn())
|
|
#Now the obvious lift of omega0 = P dx + Q dy to de Rham-Witt is [P] d[x] + [Q] d[y]
|
|
omega8_regular = regular_form(second_patch(omega8)) # The same for omega8.
|
|
omega8_regular = (second_patch(omega8_regular[0]), second_patch(omega8_regular[1]))
|
|
u = (C.x)^(-1)
|
|
v = (C.y)/(C.x)^(g+1)
|
|
omega8_lift = omega8_regular[0].teichmuller()*(u.teichmuller().diffn()) + omega8_regular[1].teichmuller()*(v.teichmuller().diffn())
|
|
aux = omega0_lift - omega8_lift - fct.teichmuller().diffn() # now aux is of the form (V(smth) + dV(smth), V(smth))
|
|
if aux.h1.function != 0:
|
|
raise ValueError('Something went wrong - aux is not of the form (V(smth) + dV(smth), V(smth)).')
|
|
decom_aux_h2 = decomposition_g0_g8(aux.h2, prec=prec) #decompose dV(smth) in aux as smth regular on U0 - smth regular on U8.
|
|
aux_h2 = decom_aux_h2[0]
|
|
aux_f = decom_aux_h2[2]
|
|
aux_omega0 = decomposition_omega0_omega8(aux.omega, prec=prec)[0]
|
|
result = superelliptic_drw_cech(omega0_lift - aux_h2.verschiebung().diffn() - aux_omega0.verschiebung(), fct.teichmuller() + aux_f.verschiebung())
|
|
return result.reduce()
|
|
|
|
def crystalline_cohomology_basis(self, prec = 50):
|
|
result = []
|
|
for a in self.de_rham_basis():
|
|
result += [de_rham_witt_lift(a, prec = prec)]
|
|
return result
|
|
|
|
superelliptic.crystalline_cohomology_basis = crystalline_cohomology_basis
|
|
|
|
def autom(self):
|
|
C = self.curve
|
|
F = C.base_ring
|
|
Rxy.<x, y> = PolynomialRing(F, 2)
|
|
Fxy = FractionField(Rxy)
|
|
if isinstance(self, superelliptic_function):
|
|
result = superelliptic_function(C, Fxy(self.function).subs({x:x+1, y:y}))
|
|
return result
|
|
if isinstance(self, superelliptic_form):
|
|
result = superelliptic_form(C, Fxy(self.form).subs({x:x+1, y:y}))
|
|
return result
|
|
if isinstance(self, superelliptic_witt):
|
|
result = superelliptic_witt(autom(self.t), autom(self.f))
|
|
return result
|
|
if isinstance(self, superelliptic_drw_form):
|
|
result = superelliptic_drw_form(0*C.x, autom(self.omega), autom(self.h2))
|
|
result += autom(self.h1).teichmuller()*(C.x + C.one).teichmuller().diffn()
|
|
return result
|
|
if isinstance(self, superelliptic_drw_cech):
|
|
result = superelliptic_drw_cech(autom(self.omega0), autom(self.f))
|
|
return result
|
|
|
|
|
|
def dy_w(C):
|
|
'''Return d[y].'''
|
|
fC = C.polynomial
|
|
fC = superelliptic_function(C, fC)
|
|
fC = fC.teichmuller()
|
|
dy_w = 1/2* ((C.y)^(-1)).teichmuller()*auxilliary_derivative(fC)
|
|
return dy_w
|
|
superelliptic.dy_w = dy_w |