A<x,y> := PolynomialRing(GF(3),2);
P,homogenize := Homogenization(A);

print("f := an improved polynomial for X1(38)");
f := x*(x-1)^3*y^14+(6*x^3-9*x^2+5*x+1)*(x-1)^2*y^13+(x-1)*(15*x^5-47*x^4+72*x^3-45*x^2+2*x+6)*y^12+(20*x^7+266*x^3-39*x^2-394*x^4-102*x^6+274*x^5-39*x+15)*y^11+(x-1)*(15*x^7-88*x^6+304*x^5-489*x^4+304*x^3+4*x^2-67*x+20)*y^10+(6*x^7-54*x^6+248*x^5-472*x^4+254*x^3+64*x^2-67*x+15)*(x-1)^2*y^9+(x-1)*(x^9-26*x^8+190*x^7-682*x^6+1049*x^5-565*x^4-117*x^3+193*x^2-50*x+6)*y^8-(x-1)*(7*x^9-73*x^8+426*x^7-981*x^6+837*x^5+14*x^4-351*x^3+126*x^2-5*x+1)*y^7-x*(x^8-14*x^7+164*x^6-444*x^5+299*x^4+240*x^3-273*x^2+22*x+17)*(x-1)^2*y^6+x*(x^8-50*x^7+197*x^6-207*x^5-177*x^4+378*x^3-107*x^2-47*x+15)*(x-1)^2*y^5-x*(x-1)*(10*x^9-61*x^8+125*x^7+24*x^6-410*x^5+490*x^4-133*x^3-91*x^2+53*x-6)*y^4-x*(x^9-6*x^8+11*x^7+40*x^6-153*x^5+139*x^4+17*x^3-63*x^2+19*x-1)*(x-1)^2*y^3-x^2*(10*x^5-32*x^4+21*x^3+25*x^2-21*x+3)*(x-1)^3*y^2-x^3*(x^4-3*x^3+x^2+8*x-3)*(x-1)^4*y-x^4*(x-1)^5 ;

f := homogenize(f);
XF2 := Scheme(ProjectiveSpace(P),f);
FF := FunctionField(XF2);
AFF := AlgorithmicFunctionField(FF);
plc1 := Places(AFF,1);
divgrp := DivisorGroup(AFF);
div1 := [divgrp ! i : i in plc1];
cuspsum := &+ div1;
load "function_degrees.m";

print("The number of places of degrees 1-6:");
[#Places(AFF,i) : i in [1..6]];

// Our goal is to show that X1(38) has no non-constant functions of degree <12
// over F3 and hence over Q. 

// Say f is a function of degree <= 11, since f^(-1)({0,1,2,infinity}) contains
// all 18 rational points, by the pigeon hole principle we may assume that
// f^(-1)({infinity} contains at least 5 distinct rational points, and so
// the remaining polar-part of f has degree at most 6.

// C0  C1 C2  C3  C4 C5  C6  C7 C8  C9 C10 C11 C12 C13 C14 C15 C16 C17 C18 C19
// [0, -3, 0, -2, -3, 0, -2, -1, 3, -2,  1,  5,  0, 0,  -1, 0,  0,  3,  2,  0]

pl,val := Support(Divisor(AFF.1));
for i in [1..#val] do
	if val[i] eq 5 then
		orb1 := pl[i];
	end if;
	if val[i] eq 2 then
		orb2 := pl[i];
	end if;
end for;

// Has one member of each orbit.
orbs := [orb1, orb2];



// Degree >1 places: [3, 3]
// Rat-places: [] (only multiplicities >1 are listed).
s := {(cuspsum + P1+P2) : P1 in Places(AFF,3), P2 in Places(AFF,3)};
print "The number of RR-computations to be done for this case is:",
#s;
s := &cat[FunctionDegrees(D) : D in s];
if not s eq [] then s := Min(s); end if;
s;

// Degree >1 places: [6]
// Rat-places: [] (only multiplicities >1 are listed).
s := {(cuspsum + P1) : P1 in Places(AFF,6)};
print "The number of RR-computations to be done for this case is:",
#s;
s := &cat[FunctionDegrees(D) : D in s];
if not s eq [] then s := Min(s); end if;
s;

// Degree >1 places: [5]
// Rat-places: [2] (only multiplicities >1 are listed).
s := {(cuspsum + P1 + 1*Q1) : P1 in Places(AFF,5), Q1 in orbs};
print "The number of RR-computations to be done for this case is:",
#s;
s := &cat[FunctionDegrees(D) : D in s];
if not s eq [] then s := Min(s); end if;
s;

// Degree >1 places: []
// Rat-places: [3, 3, 5] (only multiplicities >1 are listed).
s := {(cuspsum + 2*Q1+2*Q2+4*Q3) : Q1 in plc1, Q2 in plc1, Q3 in orbs};
print "The number of RR-computations to be done for this case is:",
#s;
s := &cat[FunctionDegrees(D) : D in s];
if not s eq [] then s := Min(s); end if;
s;

// Degree >1 places: [3]
// Rat-places: [2, 2, 2] (only multiplicities >1 are listed).
s := {(cuspsum + P1 + 1*Q1+1*Q2+1*Q3) : P1 in Places(AFF,3), Q1 in plc1, Q2 in plc1, Q3 in orbs};
print "The number of RR-computations to be done for this case is:",
#s;
s := &cat[FunctionDegrees(D) : D in s];
if not s eq [] then s := Min(s); end if;
s;

// Degree >1 places: []
// Rat-places: [2, 2, 2, 2, 3] (only multiplicities >1 are listed).
s := {(cuspsum + 1*Q1+1*Q2+1*Q3+1*Q4+2*Q5) : Q1 in plc1, Q2 in plc1, Q3 in plc1, Q4 in plc1, Q5 in orbs};
print "The number of RR-computations to be done for this case is:",
#s;
s := &cat[FunctionDegrees(D) : D in s];
if not s eq [] then s := Min(s); end if;
s;
