/////////////////////////////////////////////////////////////////////////// // Magma code that verifies the computational claims in // // M. Derickx, S. Kamienny, W. Stein, M. Stoll: // // Torsion points on elliptic curves over number fields of small degree // /////////////////////////////////////////////////////////////////////////// t0 := Cputime(); load "X1_p.magma"; // load the curves X_1(p) over F_2 P := PolynomialRing(Rationals()); // Claim after Prop. 1.3: Verify that 37 is in S(6). printf "\n=============================================================================\n"; printf "We show that 37 is in S(6).\n"; printf "=============================================================================\n\n"; E := EllipticCurve([1, 1, 1, -8, 6]); // this is 1225.b2 in the LMFDB or 1225h1 in Cremona // Set up K = Q(sqrt{5}, \cos(2 \pi/7)) K := AbsoluteField(ext); divpol37 := DivisionPolynomial(E, 37); // Factor the division polynomial fdp37 := Factorization(divpol37); assert Degree(fdp37[1,1]) eq 6; // it has factors of degree 6 // Check that the factor of degree 6 has roots in K assert not IsEmpty(Roots(fdp37[1,1], K)); // Let xi be a root. Then there is a quadratic twist E' of E over K // such that xi is the x-coordinate of a K-point of E'; this point // has order 37. printf "done.\n\n"; //========================================================================= // Theorem 3.2 // Some computations with modular symbols for X_1(29) // based on [Parent 2000]; see also William Stein's book [Stein 207], Chapter 8. p := 29; F := GF(p); p2 := ExactQuotient(p-1, 2); printf "\n=============================================================================\n"; printf "We show that J_1(%o)(Q)_tors is generated by differences of rational cusps.\n", p; printf "=============================================================================\n\n"; printf "Setting up modular symbols...\n\n"; RG := [ : b in [1..p2]] cat [ : b in F, a in [1..p2]]; ZRG := FreeAbelianGroup(#RG); function index(pair) // index in RG a := Integers()!pair[1]; if a eq 0 then b := Integers()!pair[2]; if b gt p2 then pair := ; end if; else if a gt p2 then pair := <-pair[1], -pair[2]>; end if; end if; return Position(RG, pair); end function; sigma_perm := [index(<-pair[2], pair[1]>) : pair in RG]; tau_perm := [index() : pair in RG]; sigma_ZRG := hom ZRG | [ZRG.sigma_perm[i] - ZRG.i : i in [1..#RG]]>; tau_ZRG := hom ZRG | [ZRG.tau_perm[i] - ZRG.i : i in [1..#RG]]>; fix_sigma := Kernel(sigma_ZRG); fix_tau := Kernel(tau_ZRG); // Set up the boundary map Z[R_G] --> Z[cusps] Zcusps := FreeAbelianGroup(p-1); // ordering of the cusps: (1/1, 1/2, ..., 1/p2, 1/p, 2/p, ..., p2/p) // note that 1/1 ~ 0 and 1/p ~ oo norm := func; // F mod +-1 cusp1 := func; cusp2 := func)>; // delta is the boundary map. delta := hom Zcusps | [Zcusps.cusp1(pair) - Zcusps.cusp2(pair) : pair in RG]>; ker_delta := Kernel(delta); // The relative integral homology H_1(X_1(p)(C), cusps, Z): H1_cusp, torelhom := quo; // The integral homology H_1(X_1(p)(X), Z): H1_X1p_Z, tointhom := quo; // The inclusion map of the latter into the former: incl := hom H1_cusp | [torelhom(g @@ tointhom) : g in OrderedGenerators(H1_X1p_Z)]>; printf "The relative homology H_1(X_1(%o)(C), cusps, Z) has rank %o\n", p, #Invariants(H1_cusp); printf "The homology H_1(X_1(%o)(C), Z) has rank %o\n\n", p, #Invariants(H1_X1p_Z); // Set up diamond, Hecke, star operators. // Right action on the Manin symbols of an integral matrix that is invertible mod p: act := func)]>; // Get the induced operator on the homology from an operator on Z[R_G]. induced_op_relhom := func H1_cusp | [torelhom(op(g @@ torelhom)) : g in OrderedGenerators(H1_cusp)]>>; induced_op_inthom := func H1_X1p_Z | [tointhom(op(g @@ tointhom)) : g in OrderedGenerators(H1_X1p_Z)]>>; // Extract the matrix of an endomorphism on the homology. matrix_inthom := func; matrix_relhom := func; // The diamond operator on Z[R_G]. diamond_ZRG := func ZRG | [ZRG.index() : pair in RG]>>; // The Hecke operator T_q on Z[R_G] for prime index q /= p. hecke_ZRG := func ZRG | [&+[ZRG.index(act(pair, m)) : m in mats] : pair in RG]> where mats := [Matrix(2, 2, m) : m in HeilbronnMerel(q)]>; // The star involution. // Note that this corresponds to complex conjugation on X_1(p)(C); // it is the negative of the version in Stein's book. star_ZRG := hom ZRG | [ZRG.index(<-pair[1], pair[2]>) : pair in RG]>; printf "Setting up Hecke and diamond operators...\n\n"; // Set up some of the operators on Z[R_G]. primes := PrimesInInterval(3, 7); time hecke_ops_ZRG := [hecke_ZRG(q) : q in primes]; // ~ 10 sec diamond_ops_ZRG := [diamond_ZRG(q) : q in primes]; // Determine the Eisenstein sublattice of the relative homology. // If T is a Hecke operator on the relative homology and its restriction // to the integral homology has characteristic polynomial f, then f(T) // maps the relative homology into the Eisenstein sublattice. // We know that this sublattice has rank #cusps - 1 = p - 2, // so it suffices to find one T = T_q such that the image of T has the // correct rank. In our case, T = T_3 will do. T3_relhom := induced_op_relhom(hecke_ops_ZRG[1]); T3_relhom_mat := matrix_relhom(T3_relhom); T3_inthom := induced_op_inthom(hecke_ops_ZRG[1]); cpT3_inthom := CharacteristicPolynomial(matrix_inthom(T3_inthom)); op_Eis := Evaluate(cpT3_inthom, T3_relhom_mat); assert Rank(op_Eis) eq p-2; // check that rank is correct // We get the Eisenstein sublattice by saturating the image of op_Eis. LEis := PureLattice(Lattice(op_Eis)); // We need the basis we use for the integral homology inside // the relative homology. basS := [Vector(Eltseq(incl(g))) : g in OrderedGenerators(H1_X1p_Z)]; // Now we determine the image of the integral relative homology // inside the rational homology. Its quotient by the integral homology // is (isomorphic to) the cuspidal subgroup of J_1(p)_tors. // We get the matrix that maps the integral relative homology // into the rational homology by stacking the basis of the integral // homology onto the basis of the Eisenstein sublattice, inverting // the resulting matrix and removing the last p-2 columns. basESmat := Matrix(basS cat Basis(LEis)); basESmati := ChangeRing(basESmat, Rationals())^-1; projmat := Submatrix(basESmati, 1, 1, Nrows(basESmati), #basS); // The image of the relative homology in the integral homology // (as a lattice in Q^(2 g(X_1(p)))): L0 := Lattice(IdentityMatrix(Integers(), #basS)); // integral lattice = integral homology Lcusp := Lattice(projmat) + L0; // The cuspidal subgroup of J_1(p): Qcusp, toQcusp := quo; printf "\nThe cuspidal torsion subgroup of J_1(%o)\n", p; printf "has invariants %o\n\n", Invariants(Qcusp); // Get the induced operators on the integral homology ... hecke_ops := [induced_op_inthom(T) : T in hecke_ops_ZRG]; diamond_ops := [induced_op_inthom(d) : d in diamond_ops_ZRG]; star_op := induced_op_inthom(star_ZRG); // ... and their matrices. hecke_mats := [matrix_inthom(T) : T in hecke_ops]; diamond_mats := [matrix_inthom(d) : d in diamond_ops]; star_mat := matrix_inthom(star_op); star_mat_Q := ChangeRing(star_mat, Rationals()); id_mat := Parent(star_mat)!1; // We know that for 2 < q /= p prime, the operator // T_q - q - // annihilates rational torsion; the same is true of star - 1. // We use this to find an upper bound on the rational torsion subgroup. ann_ops := [hecke_mats[j] - diamond_mats[j] - primes[j]*id_mat : j in [1..#primes]]; ann_ops_Q := [ChangeRing(mat, Rationals()) : mat in ann_ops]; L1 := &meet[Lattice(m^-1) : m in ann_ops_Q]; // Get the action of star on the quotient L1/L0 and pull back fixed points. Q1, toQ1 := quo; star_minus_1_on_Q1 := hom Q1 | [toQ1((g @@ toQ1)*star_mat_Q) - g : g in OrderedGenerators(Q1)]>; Q1fix := Kernel(star_minus_1_on_Q1); L := L0 + Lattice(Matrix([g @@ toQ1 : g in Generators(Q1fix)])); // Now L/L0 is our upper bound on J_1(p)(Q)_tors. Q, toQ := quo; printf "The group bounding the rational torsion subgroup of J_1(%o)\n", p; printf "has invariants %o\n\n", Invariants(Q); // We check that L/L0 is contained in the cuspidal torsion subgroup: assert L subset Lcusp; printf "This group is contained in the cuspidal subgroup.\n\n"; // It remains to verify that the rational (i.e., Galois invariant) part // of the cuspidal group is generated by differences of rational cusps. // The rational cusps are the cusps j/p (that map to infinity on X_0(p)). // The other cusps 1/j are permuted cyclically by Galois. // Since the Galois action commutes with the diamond operators, // the action is the same as that of the diamond operators; // in particular, it is generated by the action of <2> (for p = 29). assert Order(GF(p)!2) eq p-1; // 2 is a primitive root // We take the rational cusp oo (number p2+1) as base-point. // Images of differences cusp - oo in the rational homology: cusp_diffs := [Vector(Rationals(), Eltseq(torelhom((Zcusps.j - Zcusps.(p2+1)) @@ delta)))*projmat : j in [1..p-1]]; Zcusps_to_Qcusp := hom Qcusp | [toQcusp(Lcusp!cd) : cd in cusp_diffs]>; // The Galois action on the cusps as a permutation. // The first p2 cusps are permuted cyclically, the last p2 are fixed. Gal_action := [norm(2*F!j) : j in [1..p2]] cat [p2+j : j in [1..p2]]; Gal_on_Zcusps := hom Zcusps | [Zcusps.j : j in Gal_action]>; // The induced action on the cuspidal subgroup, minus the identity ... Gal_minus_1_on_Qcusp := hom Qcusp | [Zcusps_to_Qcusp(Gal_on_Zcusps(g @@ Zcusps_to_Qcusp)) - g : g in OrderedGenerators(Qcusp)]>; // ... and the fixed subgroup. Qcusp_fix := Kernel(Gal_minus_1_on_Qcusp); printf "The rational cuspidal subgroup has invariants %o\n", Invariants(Qcusp_fix); // Now the subgroup coming from rational cusps. Qcusp_rat := sub; printf "The subgroup coming from rational cusps has invariants %o\n\n", Invariants(Qcusp_rat); // Check they are equal. assert Qcusp_fix eq Qcusp_rat; printf "Both groups are equal!\n\n"; //========================================================================= // Proposition 3.4 and Lemma 3.7 printf "\n=============================================================================\n"; printf "We show that J_1(p)(Q) --> J_1(p)(F_2) is injective for p = 17, 29, 31, 41\n"; printf "and that assumption (b) is satisfied for these primes and d <= 7.\n"; printf "=============================================================================\n\n"; Aff2 := AffineSpace(GF(2), 2); // p = 17 printf "p = 17\n"; printf "------\n\n"; p := 17; X := X1p_models[p]; assert Genus(X) eq 5; // check that model reduces well plc := Places(X, 1); assert #plc eq ExactQuotient(p-1, 2); // sanity check printf "setting up Picard group...\n"; Cl, CltoDiv, DivtoCl := ClassGroup(X); bpt := DivtoCl(Divisor(plc[1])); printf "verifying claim...\n"; G := sub; assert Valuation(#G, 2) eq 3; // see Table 1 in [CES03] printf "successful!\n\n"; // p = 29 printf "p = 29\n"; printf "------\n\n"; p := 29; X := X1p_models[p]; assert Genus(X) eq 22; // check that model reduces well plc := Places(X, 1); assert #plc eq ExactQuotient(p-1, 2); // sanity check printf "setting up Picard group...\n"; time Cl, CltoDiv, DivtoCl := ClassGroup(X); // takes a bit more than one minute bpt := DivtoCl(Divisor(plc[1])); printf "verifying first claim...\n"; G := sub; assert Valuation(#G, 2) eq 6; // see Table 1 in [CES03]; note that the entry there must be divided by 2^6 // Check assumption (b) // no points of degree <= 5: printf "verifying second claim...\n"; assert forall{d : d in [2,3,4,5] | IsEmpty(Places(X, d))}; plc6 := Places(X, 6); plc7 := Places(X, 7); // Check that the prime divisors of degree 6 or 7 do not map into G. assert forall{pl : pl in plc6 | DivtoCl(Divisor(pl)) - 6*bpt notin G}; assert forall{pl : pl in plc7 | DivtoCl(Divisor(pl)) - 7*bpt notin G}; printf "successful!\n\n"; // p = 31 printf "p = 31\n"; printf "------\n\n"; p := 31; X := X1p_models[p]; assert Genus(X) eq 26; // check that model reduces well time plc := Places(X, 1); assert #plc eq ExactQuotient(p-1, 2); // sanity check printf "setting up Picard group...\n"; time Cl, CltoDiv, DivtoCl := ClassGroup(X); // takes about 1.5 minutes bpt := DivtoCl(Divisor(plc[1])); printf "verifying first claim...\n"; G := sub; assert Valuation(#G, 2) eq 2; // see Table 1 in [CES03] // Check assumption (b) // no points of degree <= 4: printf "verifying second claim...\n"; assert forall{d : d in [2,3,4] | IsEmpty(Places(X, d))}; plc5 := Places(X, 5); plc6 := Places(X, 6); plc7 := Places(X, 7); // Check that the prime divisors of degrees 5 to 7 do not map into G. assert forall{pl : pl in plc5 | DivtoCl(Divisor(pl)) - 5*bpt notin G}; assert forall{pl : pl in plc6 | DivtoCl(Divisor(pl)) - 6*bpt notin G}; assert forall{pl : pl in plc7 | DivtoCl(Divisor(pl)) - 7*bpt notin G}; printf "successful!\n\n"; // p = 41 printf "p = 41\n"; printf "------\n\n"; p := 41; X := X1p_models[p]; assert Genus(X) eq 51; // check that model reduces well time plc := Places(X, 1); // less than half a minute assert #plc eq ExactQuotient(p-1, 2); // sanity check printf "setting up Picard group (this will take several hours)...\n"; time Cl, CltoDiv, DivtoCl := ClassGroup(X); // takes about 8 hours bpt := DivtoCl(Divisor(plc[1])); printf "verifying first claim...\n"; G := sub; assert Valuation(#G, 2) eq 4; // see Table 1 in [CES03] // Check assumption (b) printf "verifying second claim...\n"; // no points of degree <= 7 except in degree 5: assert forall{d : d in [2,3,4,6,7] | IsEmpty(Places(X, d))}; // but two orbits of points of degree 5: plc5 := Places(X, 5); assert #plc5 eq 2; // Check that the prime divisors of degree 5 do not map into G. assert forall{pl : pl in plc5 | DivtoCl(Divisor(pl)) - 5*bpt notin G}; printf "successful!\n\n"; //========================================================================= // Lemma 5.4 printf "\n=============================================================================\n"; printf "We prove Lemma 5.4.\n"; printf "=============================================================================\n\n"; // Kamienny's criterion for X_0(p) function kamienny0(p, l, d : max_n := 60, max_q := 20) // Set up modular symbols for Gamma_0(p) M := ModularSymbols(Gamma0(p)); MS := CuspidalSubspace(M); // the winding element as an element of the rational homology e := WindingElementProjection(MS); S := ModularSymbols(p, 2, -1); // dual to the space of cusp forms function t1(n) // apply t_1 to t_0 = T_n Tn_MS := HeckeOperator(MS, n); Tn_S := IntegralHeckeOperator(S, n); pol := CharacteristicPolynomial(Tn_S); fpol := Factorization(pol); factors := [f : f in fpol | f[2] ge 2 or e*Evaluate(ExactQuotient(pol, f[1]^f[2]), Tn_MS) eq 0]; // printf "exponents: %o\n", [f[2] : f in fpol]; return Evaluate(&*[Parent(pol) | f[1]^f[2] : f in factors], Tn_S); end function; id := IntegralHeckeOperator(S, 1); if l eq 2 and p mod 8 eq 1 then t2_list := [ : q in PrimesInInterval(3, max_q) | q ne p]; else t2_list := [<1, id>]; end if; // Hecke operators T_1, ..., T_d mod l hecke := [ChangeRing(IntegralHeckeOperator(S, j), GF(l)) : j in [1..d]]; for n := 2 to max_n do t1n := t1(n); for pair in t2_list do t := ChangeRing(t1n*pair[2], GF(l)); // check linear independence of t*T_j for j = 1,...,d if Rank(Matrix([Eltseq(t*T) : T in hecke])) eq d then return true, n, pair[1]; end if; end for; end for; return false, _, _; end function; function kamienny0rc(p, l, d : max_n := 60, max_q := 20, ops := [2,3,5,7]) // Set up modular symbols for Gamma_0(p) M := ModularSymbols(Gamma0(p)); MS := CuspidalSubspace(M); // the winding element as an element of the rational homology e := WindingElementProjection(MS); S := ModularSymbols(p, 2, -1); // dual to the space of cusp forms Tnseq_S := [IntegralHeckeOperator(S, n) : n in ops]; Tnseq_MS := [HeckeOperator(MS, n) : n in ops]; function t1(cofs) // apply t_1 to t_0 = cofs[1]*T_2 + ... + cofs[4]*T_7 Tn_MS := &+[cofs[i]*Tnseq_MS[i] : i in [1..#ops]]; Tn_S := &+[cofs[i]*Tnseq_S[i] : i in [1..#ops]]; pol := CharacteristicPolynomial(Tn_S); fpol := Factorization(pol); factors := [f : f in fpol | f[2] ge 2 or e*Evaluate(ExactQuotient(pol, f[1]^f[2]), Tn_MS) eq 0]; // printf "exponents: %o\n", [f[2] : f in fpol]; return Evaluate(&*[Parent(pol) | f[1]^f[2] : f in factors], Tn_S); end function; id := IntegralHeckeOperator(S, 1); if l eq 2 and p mod 8 eq 1 then t2_list := [ : q in PrimesInInterval(3, max_q) | q ne p]; else t2_list := [<1, id>]; end if; // Hecke operators T_1, ..., T_d mod l hecke := [ChangeRing(IntegralHeckeOperator(S, j), GF(l)) : j in [1..d]]; for n := 1 to max_n do cofs := [Random({-1..1}) : i in [1..#ops]]; t1n := t1(cofs); for pair in t2_list do t := ChangeRing(t1n*pair[2], GF(l)); // check linear independence of t*T_j for j = 1,...,d if Rank(Matrix([Eltseq(t*T) : T in hecke])) eq d then return true, cofs, pair[1]; end if; end for; end for; return false, _, _; end function; // Run the check for all primes in plist. // Returns the sequence of primes, for which the criterion failed // and the lagrest n and q used for those where it was successful. function check_kam0(plist, l, d) max_n := 0; max_q := 0; fail := []; for p in plist do t0 := Cputime(); flag, n, q := kamienny0(p, l, d); t := Cputime(t0); if flag then max_n := Max(max_n, n); max_q := Max(max_q, q); if l eq 2 then printf "p = %o: success with T_%o and q = %o (%o sec)\n", p, n, q, t; else printf "p = %o: success with T_%o (%o sec)\n", p, n, t; end if; else Append(~fail, p); printf "p = %o: failure (%o sec)\n", p, t; end if; end for; return fail, max_n, max_q; end function; // Check the criterion time fail7, maxn7, maxq7 := check_kam0(PrimesInInterval(17, 2281), 2, 7); // Time: 5367.670 printf "\nd = 7 --> max. n used: %o, max. q used: %o\n\n", maxn7, maxq7; // 14, 5 assert fail7 eq PrimesInInterval(17, 127) cat [137, 157, 193, 431]; // use random linear combination of T_2,T_3,T_5,T_7 as t_0 for p = 431 assert kamienny0rc(431, 2, 7); // T_2 + T_3 - T_7 works. Prune(~fail7); // remove 431 from fail7 time fail6, maxn6, maxq6 := check_kam0(fail7, 2, 6); // Time: 1.170 printf "\nd = 6 --> max. n used: %o, max. q used: %o\n\n", maxn6, maxq6; // 3, 3 assert fail6 eq PrimesInInterval(17, 67) cat [73] cat PrimesInInterval(79, 103) cat [109, 113, 193]; time fail5, maxn5, maxq5 := check_kam0(fail6, 2, 5); // Time: 0.630 printf "\nd = 5 --> max. n used: %o, max. q used: %o\n\n", maxn5, maxq5; // 10, 3 assert fail5 eq PrimesInInterval(17, 53) cat [61, 67, 73, 79, 89, 97, 101]; time fail4, maxn4, maxq4 := check_kam0(fail5, 2, 4); // Time: 0.440 printf "\nd = 4 --> max. n used: %o, max. q used: %o\n\n", maxn4, maxq4; // 3, 3 assert fail4 eq PrimesInInterval(17, 43) cat [53, 61, 67, 73, 79, 97, 101]; time fail3, maxn3, maxq3 := check_kam0(fail4, 2, 3); // Time: 0.300 printf "\nd = 3 --> max. n used: %o, max. q used: %o\n\n", maxn3, maxq3; // 5, 3 assert fail3 eq PrimesInInterval(17, 43) cat [73, 79]; printf "\nsuccessful!\n\n"; printf "\n=============================================================================\n"; printf "We prove Lemma 5.6.\n"; printf "=============================================================================\n\n"; // Version for X_H with speed-up. function kamiennyH(p, d, h : max_n := 60, max_q := 20) // h represents a generator of H. // Set up modular symbols for Gamma_H. M := ModularSymbolsH(p, [h], 2, 0); MS := CuspidalSubspace(M); e := WindingElementProjection(MS); S := ModularSymbolsH(p, [h], 2, -1); // dual to the space of cusp forms function diamond(q) // diamond operator for q prime op := IntegralHeckeOperator(S, q)^2 - IntegralHeckeOperator(S, q^2); return ChangeRing(1/q*ChangeRing(op, Rationals()), Integers()); end function; function t1(n) // apply t_1 to t_0 = T_n Tn_MS := HeckeOperator(MS, n); Tn_S := IntegralHeckeOperator(S, n); pol := CharacteristicPolynomial(Tn_S); fpol := Factorization(pol); factors := [f : f in fpol | f[2] ge 2 or e*Evaluate(ExactQuotient(pol, f[1]^f[2]), Tn_MS) eq 0]; // printf "exponents: %o\n", [f[2] : f in fpol]; return Evaluate(&*[Parent(pol) | f[1]^f[2] : f in factors], Tn_S); end function; id := IntegralHeckeOperator(S, 1); // We could try to figure out whether the rational torsion subgroup // on the winding quotient has odd order, in which case // we can set t2_list := [<1, id>]; t2_list := [ : q in PrimesInInterval(3, max_q) | q ne p]; // Hecke operators T_1, ..., T_d mod 2 hecke := [ChangeRing(IntegralHeckeOperator(S, j), GF(2)) : j in [1..d]]; r := 1; repeat r := NextPrime(r); until Order(GF(p)!r) eq p-1; // r is prime and generates the quotient (Z/pZ)^*/+-H dia := ChangeRing(diamond(r), GF(2)); U, fromU := UnitGroup(GF(p)); index := Index(U, sub); dias := [dia^j : j in [0..index-1]]; d2 := d div 2; for n := 2 to max_n do vprintf User1: "Working with t_0 = T_%o\n", n; t1n := t1(n); for pair in t2_list do vprintf User1: " Using second factor with q = %o\n", pair[1]; t := ChangeRing(t1n*pair[2], GF(2)); products := [[tt*T : T in hecke[1..j eq 1 select d else d2]] where tt := t*dias[j] : j in [1..#dias]]; mat := Matrix([Eltseq(op) : op in &cat products]); ker := Kernel(mat); vprintf User1: " dim(ker) = %o\n", Dimension(ker); // Use linear codes machinery to enumerate short relations. code := LinearCode(ker); dist := PartialWeightDistribution(code, d); num := &+[Integers() | e[2] : e in dist[2..#dist]]; // number of nonzero words // Skip if number is too large if num gt 10000 then vprintf User1: " too many short relations\n"; continue pair; // try next pair of operators else rels := &join{Words(code, l) : l in [1..d]}; end if; for rel in rels do cofs := Eltseq(rel); // Check if relation is of the correct (i.e., forbidden) form parts := [cofs[1..d]] cat [cofs[d+j*d2+1..d+(j+1)*d2] : j in [0..index-2]]; // The sum of the (1-based) position numbers of the last "1" // in each part must be <= d. positions := [last_one(part) : part in parts] where last_one := function(seq) for i := #seq to 1 by -1 do if seq[i] eq GF(2)!1 then return i; end if; end for; return 0; end function; if &+positions le d then // offending relation found vprintf User1: " linear depenence found!\n"; continue pair; // try next pair of operators end if; end for; // If we get here, no offending relations were found ==> success. return true, n, pair[1]; end for; end for; return false, _, _; end function; // Determine a representative of a generator of H from H's index. function genH(p, index) assert IsDivisibleBy(p-1, 2*index); r := PrimitiveRoot(p); return Integers()!((GF(p)!r)^index); end function; function kamienny(p, d) // Try all possible H's, starting with smallest index > 1 // (we assume we already did kamienny0). t0 := Cputime(); indices := Divisors(ExactQuotient(p-1, 2)); for ind in indices[2..#indices] do flag, n, q := kamiennyH(p, d, genH(p, ind)); if flag then t := Cputime(t0); printf "p = %o: success with index(H) = %o, T_%o and q = %o (%o sec)\n", p, ind, n, q, t; return true, ind, n, q; end if; end for; return false, _, _, _; end function; // Test the criterion for the primes left over from before. assert forall{p : p in fail7 | p lt 37 or kamienny(p, 7)}; assert forall{p : p in [23, 31, 37] | kamienny(p, 6)}; assert kamienny(19, 4); assert kamienny(29, 3); // A separate check that we use in the proof of Lemma 7.3 (127 notin S(7)). assert kamiennyH(127, 7, 28); printf "\nsuccessful!\n\n"; //========================================================================= printf "\n=============================================================================\n"; printf "We prove Theorem 6.1, based on Corollary 6.8.\n"; printf "=============================================================================\n\n"; // The list of primes M that can be used to verify the criterion of Cor. 6.8. // for d = (1, 2,) 3, 4, ..., 25. Mseq := [0, 0, 29, 37, 41, 43, 47, 47, 53, 53, 53, 61, 73, 73, 79, 79, 89, 89, 89, 101, 101, 109, 109, 109, 127]; // The function epsilon_M(a). function eps(M, a) r := Integers()!a; return 2*r lt M select 0 else 1; end function; // This checks the criterion for a given d with M = Mseq[d]. function check_oesterle(d) M := Mseq[d]; for u in GF(M) do mat := Matrix(GF(M), [[eps(M,n*a) - eps(M,n*u/a) : a in GF(M) | a ne 0] : n in [1..d]]); if Rank(mat) lt d then return false; end if; end for; return true; end function; for d := 3 to 25 do assert check_oesterle(d); end for; printf "successful!\n\n"; //========================================================================= printf "\n=============================================================================\n"; printf "We prove Corollary 7.2.\n"; printf "=============================================================================\n\n"; // We need to check for p = 113, 127 that the positive-rank simple factors // of J_1(p) all occur in J_0(p). // We can check whether the simple factor associated to a Galois orbit // of newforms has positive rank (<==> L(*,1) = 0) in Magma via // IsZero(WindingElementProjection(M)) . assert #[M : M in NewformDecomposition(CuspidalSubspace(ModularSymbols(Gamma1(113),2,1))) | IsZero(WindingElementProjection(M))] eq #[M : M in NewformDecomposition(CuspidalSubspace(ModularSymbols(Gamma0(113),2,1))) | IsZero(WindingElementProjection(M))]; assert #[M : M in NewformDecomposition(CuspidalSubspace(ModularSymbols(Gamma1(127),2,1))) | IsZero(WindingElementProjection(M))] eq #[M : M in NewformDecomposition(CuspidalSubspace(ModularSymbols(Gamma0(127),2,1))) | IsZero(WindingElementProjection(M))]; // We also check the inequalities. assert 975/48/4096*(113^2-1)/7 gt 9 and 975/48/4096*(71^2-1)/3 gt 8; printf "successful!\n\n"; //========================================================================= printf "\n=============================================================================\n"; printf "We prove Lemma 8.1.\n"; printf "=============================================================================\n\n"; // Brute-force verification of Waterhouse's theorem: // This constructs representatives E of all isomorphism classes // of elliptic curves over the finite field F and returns // the orders of the groups E(F) as a set. orders := func; // primes(F) returns the set of all prime divisors of the orders. // This is the same as the set of all possible prime orders // of F-points on an elliptic curve over F. // Equivalently, it is the set of primes p such that // X_1(p)(F) contains non-cuspidal points. primes := func; prlists := [primes(GF(2, d)) : d in [1..7]]; prlists := [&join prlists[1..d] : d in [1..7]]; // Verify Lemma 7.1. // Note that 41 is not in primes(GF(2,d)) for d = 6,7, // but is excluded, since 41 occurs for d = 5. assert prlists[3] eq Set(PrimesUpTo(7)) join {13}; assert prlists[4] eq Set(PrimesUpTo(17)); assert prlists[5] eq Set(PrimesUpTo(19)) join {41}; assert prlists[6] eq (Set(PrimesUpTo(37)) diff {23}) join {41, 73}; assert prlists[7] eq Set(PrimesUpTo(43)) join {59, 61, 67, 71, 73, 113}; // 31 for d >= 5 and 127 for d >= 7 are excluded, but this is related to cusps. printf "successful!\n\n"; //========================================================================= // Lemmas 8.4 and 8.5 // Partition of point sets into orbits under Frobenius. frob := func; function frob_orbits(pts) orbits := [Parent({Rep(pts)})|]; pts := Set(pts); // pts is an indexed set; turn into a set while not IsEmpty(pts) do rep := Rep(pts); orb := {rep}; while true do rep := frob(rep); if rep in orb then break; end if; Include(~orb, rep); end while; Append(~orbits, orb); pts diff:= orb; end while; return orbits; end function; printf "\n=============================================================================\n"; printf "We prove Lemma 8.4: assumption (b) is satisfied for (d,p) = (6,73).\n"; printf "=============================================================================\n\n"; // Precision bound for q-expansions prec := 500; printf "setting up modular symbols...\n\n"; // Set up modular symbols for Gamma_H with H = <10> in (Z/73Z)^*/{+-1}. // The following function is not mentioned in the Magma documentation. // It sets up the space of modular symbols associated to Gamma_H; // see .../package/Geometry/ModSym/multichar.m in the Magma directory. M := ModularSymbolsH(73, [10], 2, 1); MS := CuspidalSubspace(M); g := Dimension(MS); // 43, the genus of X_H // Get a basis of differentials with integral q-expansions. qbas := qIntegralBasis(MS, prec); // Find quadratic relations between them. printf "constructing canonical model of X_H over F_2...\n\n"; P := PolynomialRing(Integers(), g); mons2 := MonomialsOfDegree(P, 2); monsq := [Evaluate(m, qbas) : m in mons2]; mat := Matrix([[Coefficient(s, j) : j in [2..prec]] : s in monsq]); kermat := KernelMatrix(mat); // takes a few seconds pols := [&+[kermat[i,j]*mons2[j] : j in [1..#mons2]] : i in [1..Nrows(kermat)]]; // Check that the space of quadrics vanishing on the image // has the correct dimension. // Note: binom(g+1, 2) is number of monomials of degree 2, // the dimension of L(2*canonical) is 3(g-1). assert #pols eq Binomial(g+1, 2) - 3*(g-1); // Reduce the model of X_H mod 2. PF2 := PolynomialRing(GF(2), g); pols2 := [PF2!pol : pol in pols]; Pr := ProjectiveSpace(PF2); X := Curve(Pr, pols2); printf "find points over F_{2^6}...\n\n"; // We want points over F_{2^6}, so base-change to that field. X_6 := BaseChange(X, GF(2^6)); // Find all points on X over that field. time pts6 := &join{Points(Scheme(X_6, X_6.2-a*X_6.1)) : a in GF(2^6)} join Points(Scheme(X_6, X_6.1)); // about 6 minutes // Partition into orbits under Frobenius. orbits := frob_orbits(pts6); printf "determine factor of winding quotient...\n\n"; // Extract the 30-dimensional simple factor of J_H. decomp := NewformDecomposition(MS); assert [Dimension(D) : D in decomp] eq [1,2,2,8,30]; S30 := decomp[5]; // Verify that the associated abelian variety belongs to the // winding quotient and all isogeneous abelian varieties // have odd order torsion. assert not IsZero(WindingElementProjection(S30)); assert IsOdd(TorsionBound(S30, 10)); printf "set up projection...\n\n"; // Find the corresponding projection phi: X --> P^{dim A - 1}. // Get Z-basis of differentials in terms of q-expansions. qbasS := qIntegralBasis(S30, prec); // Write them as linear combinations on the basis of MS. matMS := Matrix([[Coefficient(s, j) : j in [1..prec-1]] : s in qbas]); qbasSinqbas := Matrix([Solution(matMS, Vector([Coefficient(b, j) : j in [1..prec-1]])) : b in qbasS]); // This gives the linear forms over F_2 realizing the projection phi. linforms := [&+[qbasSinqbas[i,j]*PF2.j : j in [1..Ncols(qbasSinqbas)]] : i in [1..Nrows(qbasSinqbas)]]; // Pick representatives of the honest degree 6 points. reps := [Rep(o) : o in orbits | #o eq 6]; printf "verify criterion...\n\n"; // In principle, we only need to check the conditions on the one // orbit that is the image of an orbit of degree-6 points on X_1(73), // but it is hard to figure out which one it is. // So we do it for all representative points we found. // Evaluate at each point in reps to check that not all are zero. lfeval := [[Evaluate(l, Eltseq(pt)) : l in linforms] : pt in reps]; assert forall{a : a in lfeval | exists{b : b in a | b ne 0}}; // Check for each point that its coordinates generate a 6-dim'l space over F_2. assert forall{a : a in lfeval | Rank(Matrix([Eltseq(b) : b in a])) eq 6}; printf "successful!\n\n"; printf "verify that there is a rational point on X_H^(6) in the relevant residue class\n\n"; // Verify that there is a rational point on X_H^(6) reducing to x_H. // Set up field of degree 24. P := PolynomialRing(Rationals()); f := 73*x^24 + 222504*x^23 + 248310960*x^22 - 16333087744*x^21 + 24295278336*x^20 + 37297672387584*x^19 + 1250056453574656*x^18 + 3538765245136896*x^17 + 14951152026648576*x^16 - 44656705696694272*x^15 + 3479952968122368*x^14 + 3688842876847915008*x^13 + 20088509961884139520*x^12 + 43140842422831939584*x^11 + 48643981386046242816*x^10 + 37264018188680560640*x^9 + 31484738975758811136*x^8 + 27982614873004572672*x^7 + 19111830342294568960*x^6 + 12377884091083653120*x^5 + 6744642918801211392*x^4 + 1626555929480658944*x^3 + 83861950873731072*x^2 - 12244161486913536*x + 281474976710656; K := NumberField(f); // Define an elliptic curve over K with CM by -3. EK := EllipticCurve([0, a^2*(a+1)^3]); // Verify that it has a K-point of order 73. ptEK := EK![a*(a+1), a*(a+1)^2]; assert 73*ptEK eq EK!0; // We now show that 2 splits into four primes of degree 6 in K. OK := MaximalOrder(GeneralisedEquationOrder(K)); // takes a few seconds dec2 := Decomposition(OK, 2); assert #dec2 eq 4 and forall{e : e in dec2 | Degree(e[1]) eq 6 and e[2] eq 1}; // This implies that the 24 geometric points over F_2-bar that we // obtain as reductions of the 24 geometric points over Q-bar coming // from (EK, ptEK) on X_1(73) fall into four Galois orbits of size 6. // It follows that these points are the points mapping into the support // of x_H, since they are the only non-cuspidal points of degree 6. printf "successful!\n\n"; //========================================================================= printf "\n=============================================================================\n"; printf "We prove Lemma 8.5: assumption (b) is satisfied for (d,p) = (7,43).\n"; printf "=============================================================================\n\n"; // Check condition on alpha' for d = 7, p = 43 d := 7; p := 43; // Precision bound for q-expansions prec := 2000; printf "setting up modular symbols...\n\n"; // Set up modular symbols for Gamma_H. M := ModularSymbolsH(p, [1], 2, 1); MS := CuspidalSubspace(M); g := Dimension(MS); // the genus of X_H // Get a basis of differentials with integral q-expansions. qbas := qIntegralBasis(MS, prec); // Find quadratic relations between them. printf "constructing canonical model of X_1(%o) over F_2...\n\n", p; P := PolynomialRing(Integers(), g); mons2 := MonomialsOfDegree(P, 2); monsq := [Evaluate(m, qbas) : m in mons2]; mat := Matrix([[Coefficient(s, j) : j in [2..prec]] : s in monsq]); kermat := KernelMatrix(mat); // takes a few seconds pols := [&+[kermat[i,j]*mons2[j] : j in [1..#mons2]] : i in [1..Nrows(kermat)]]; // Check that the space of quadrics vanishing on the image // has the correct dimension. // Note: binom(g+1, 2) is number of monomials of degree 2, // the dimension of L(2*canonical) is 3(g-1). assert #pols eq Binomial(g+1, 2) - 3*(g-1); matMS := Matrix([[Coefficient(s, j) : j in [1..prec-1]] : s in qbas]); // Reduce the model of X_H mod 2. PF2 := PolynomialRing(GF(2), g); pols2 := [PF2!pol : pol in pols]; Pr := ProjectiveSpace(PF2); X := Curve(Pr, pols2); // We want points over F_{2^d}, so base-change to that field. X_d := BaseChange(X, GF(2^d)); printf "setting up modular symbols...\n\n"; // Find projection to another X_H S7 := CuspidalSubspace(ModularSymbolsH(p, [genH(p, 7)], 2, 1)); qbasS7 := qIntegralBasis(S7, prec); // Write them as linear combinations on the basis of MS. qbasS7inqbas := Matrix([Solution(matMS, Vector([Coefficient(b, j) : j in [1..prec-1]])) : b in qbasS7]); // This gives the linear forms over F_2 realizing the projection phi. linforms7 := [&+[qbasS7inqbas[i,j]*PF2.j : j in [1..Ncols(qbasS7inqbas)]] : i in [1..Nrows(qbasS7inqbas)]]; g7 := Dimension(S7); printf "constructing canonical model of X_H over F_2...\n\n"; P7 := PolynomialRing(Integers(), g7); mons2_7 := MonomialsOfDegree(P7, 2); monsqS7 := [Evaluate(m, qbasS7) : m in mons2_7]; matS7 := Matrix([[Coefficient(s, j) : j in [2..prec]] : s in monsqS7]); kermatS7 := KernelMatrix(matS7); // takes a few seconds polsS7 := [&+[kermatS7[i,j]*mons2_7[j] : j in [1..#mons2_7]] : i in [1..Nrows(kermatS7)]]; assert #polsS7 eq Binomial(g7+1, 2) - 3*(g7-1); PF2_7 := PolynomialRing(GF(2), g7); pols2S7 := [PF2_7!pol : pol in polsS7]; Pr7 := ProjectiveSpace(PF2_7); X7 := Curve(Pr7, pols2S7); X7_d := BaseChange(X7, GF(2^d)); matMS7 := Matrix([[Coefficient(s, j) : j in [1..prec-1]] : s in qbasS7]); printf "find points over F_{2^%o} on X_H...\n\n", d; // Find all points on X7 over F_{2^d}. time pts7 := &join{Points(Scheme(X7_d, X7_d.2-a*X7_d.1)) : a in GF(2^d)} join Points(Scheme(X7_d, X7_d.1)); // Partition in to orbits under Frobenius. orbits7 := frob_orbits(pts7); // The map X_1(p) --> X_H is of degree 3, prime to d, // so the points that interest us are lifts of honest degree d points. reps7 := [Rep(o) : o in orbits7 | #o eq d]; printf "reduce to points lifting to X_1(%o)...\n\n", p; // Find all points on X over F_{2^d} by lifting these points time points := [* Points(Scheme(X_d, Minors(Matrix([[Evaluate(l, [X_d.i : i in [1..g]]) : l in linforms7], Eltseq(pt)]), 2))) : pt in reps7 *]; // Representatives of he points on X_H that come from degree d points on X honest_pts7 := [reps7[i] : i in [1..#reps7] | not IsEmpty(points[i])]; pts := &join[s : s in points | not IsEmpty(s)]; printf "determine factor of winding quotient...\n\n"; // Extract simple factors of J_H. decomp := NewformDecomposition(S7); assert [Dimension(D) : D in decomp] eq [1, 2, 6, 6]; S := decomp[2] + decomp[3] + decomp[4]; // Verify that the associated abelian variety belongs to the // winding quotient and all isogeneous abelian varieties // have odd order torsion. assert forall{j : j in [2,3,4] | not IsZero(WindingElementProjection(decomp[j]))}; assert IsOdd(TorsionBound(S, 20)); printf "set up projection...\n\n"; // Find the corresponding projection phi: X --> P^{dim A - 1}. // Get Z-basis of differentials in terms of q-expansions. qbasS := qIntegralBasis(S, prec); // Write them as linear combinations on the basis of MS. qbasSinqbasS7 := Matrix([Solution(matMS7, Vector([Coefficient(b, j) : j in [1..prec-1]])) : b in qbasS]); // This gives the linear forms over F_2 realizing the projection phi. linforms := [&+[qbasSinqbasS7[i,j]*PF2_7.j : j in [1..Ncols(qbasSinqbasS7)]] : i in [1..Nrows(qbasSinqbasS7)]]; printf "verify criterion...\n\n"; // Pick representatives of the honest degree d points // (we did that before on X_H). // There are six reps, three come from elliptic curves and three from // non-rational cusps. reps := honest_pts7; // it also works with reps7, so we do not need to use X_1(p) // Evaluate at each point in reps to check that not all are zero. lfeval := [[Evaluate(l, Eltseq(pt)) : l in linforms] : pt in reps]; assert forall{a : a in lfeval | exists{b : b in a | b ne 0}}; // Check for each point that its coordinates generate a d-dim'l space over F_2. assert forall{a : a in lfeval | Rank(Matrix([Eltseq(b) : b in a])) eq d}; printf "successful!\n\n"; //========================================================================= printf "\n=============================================================================\n"; printf "We prove Lemma 8.7: assumption (b) holds for d = 7 and p = 59, 61, 67, 73.\n"; printf "=============================================================================\n\n"; // For p = 59, 67, 73, check that all positive-rank simple factors // of J_1(p) occur in J_0(p). printf "checking claims on positive-rank simple factors of J_1(p)...\n\n"; for p in [59, 67, 73] do assert #[M : M in NewformDecomposition(CuspidalSubspace(ModularSymbols(Gamma1(p),2,1))) | IsZero(WindingElementProjection(M))] eq #[M : M in NewformDecomposition(CuspidalSubspace(ModularSymbols(Gamma0(p),2,1))) | IsZero(WindingElementProjection(M))]; end for; // For p = 61, check that all positive-rank simple factors occur in J_{<3>}. assert #[M : M in NewformDecomposition(CuspidalSubspace(ModularSymbols(Gamma1(61),2,1))) | IsZero(WindingElementProjection(M))] eq #[M : M in NewformDecomposition(CuspidalSubspace(ModularSymbolsH(61,[3],2,1))) | IsZero(WindingElementProjection(M))]; printf "sucessful!\n\n"; // For the further computations, we first set up lists of all elliptic // curves (up to isomorphism) over the fields F_{2^n}, for n = 1, ..., 7, // restricting to one representative of each orbit of Forbenius. function F2n_mod_frob(n) F := GF(2, n); set := Set(F); result := {F| }; while not IsEmpty(set) do a := Rep(set); Include(~result, a); b := a; repeat Exclude(~set, b); b := b^2; until b eq a; end while; return result; end function; curves := [[ : E in &cat[Twists(EllipticCurveFromjInvariant(j)) : j in F2n_mod_frob(n)]] : n in [1..7]]; // We translate between non-cuspidal points on X_1(p) over F_2 // and pairs (E, P) consisting of an elliptic curve E over F = F_{2^n} // and a point P in E(F) of exact order p by representing (E, P) by // its Tate parameters, i.e., (b, c) in F^2 such that (E, P) is isomorphic to // (Y^2 + (1-c) XY - b Y = X^3 - b X^2, (0, 0)) . // The coordinates (x, y) of the models of X_1(p) we use // are related to (b, c) via // r = b/c, s = c^2/(b-c); b = rs(r-1), c = s(r-1) // x = (s-r)/(rs-2r+1), y = (rs-2r+1)/(s^2-s-r+1); // r = (x^2y-xy+y-1)/(x(xy-1)), s = (xy-y+1)/(xy). function EP_to_xy(E, P) // E in elliptic curve over some field F, P in E(F) not of order 1, 2, 3. // Returns (x, y) in F^2 derived from Tate parameters (b, c) as above. // First transform to Tate form. assert P[3] eq 1; // P is not the origin // extract affine coordinates of P x0 := P[1]; y0 := P[2]; // and the coefficients of E a1, a2, a3, a4, a6 := Explode(aInvariants(E)); // shift coordinates to put P at (0, 0) a1, a2, a3, a4 := Explode([a1, a2 + 3*x0, a3 + a1*x0 + 2*y0, a4 + 2*x0*a2 - a1*y0 + 3*x0^2]); // make y = 0 the tangent line at (0, 0) lambda := a4/a3; // slope of tangent a1, a2, a3 := Explode([a1 + 2*lambda, a2 - lambda*a1 - lambda^2, a3]); // scale x and y to make a2 = a3 --> get b = -a2, c = 1 - a1 u := a2/a3; // scaling factor b := -a2*u^2; c := 1 - a1*u; // compute x and y from b and c r := b/c; s := c^2/(b - c); t := r*s - 2*r + 1; x := (s - r)/t; y := t/(s^2 - s - r + 1); return x, y; end function; function xy_to_EP(x, y) // x, y in a field F. // Returns E in Tate form corresponding to (x, y) (with P = E![0,0]). // compute b and c from x and y xy := x*y; t := xy - y + 1; s := t/xy; r := (x*xy - t)/(x*(xy - 1)); c := s*(r - 1); b := r*c; // construct Tate curve return EllipticCurve([1-c, -b, -b, 0, 0]); end function; function place_from_point(pt) // pt is a point on a curve over an extension of the base field. // Returns the unique place determined by pt. plcs := Places(pt); assert #plcs eq 1; return plcs[1]; end function; // Compute the effect of a diamond operator on (x, y). function diamondop(a, x, y) E := xy_to_EP(x, y); return EP_to_xy(E, a*E![0,0]); end function; // Now the effect on a divisor. function diamondop_div(a, D) // split divisor into places with multiplicities plcs, mults := Support(D); X := Curve(D); // apply to a point of each place and turn back into places pts := [* Eltseq(RepresentativePoint(plc)) : plc in plcs *]; assert forall{1 : s in pts | s[3] eq 1}; // all affine points pts := [X(Universe(s))![x,y,1] where x, y := diamondop(a, s[1], s[2]) : s in pts]; newplcs := [place_from_point(pt) : pt in pts]; assert forall{i : i in [1..#plcs] | Degree(plcs[i]) eq Degree(newplcs[i])}; // sanity check return &+[Parent(D)| mults[i]*newplcs[i] : i in [1..#plcs]]; end function; // Given x, y, return the sequence of four pairs // corresponding to 3-isogenous curves with the image of P. function T3(x, y) E := xy_to_EP(x, y); pol := DivisionPolynomial(E, 3); // each linear factor corresponds to one isogeny F := SplittingField(pol); fact := Factorization(ChangeRing(pol, F)); lins := &cat[[e[1] : i in [1..e[2]]] : e in fact]; EF := BaseChange(E, F); isogs := [ where EE, isog := IsogenyFromKernel(EF, lin) : lin in lins]; return [ where x1, y1 := EP_to_xy(e[1], e[2](EF![0,0])) : e in isogs]; end function; // The effect of T_3 on a divisor. function T3_div(D) // split divisor into places with multiplicities plcs, mults := Support(D); X := Curve(D); // collect result result := Parent(D)!0; // apply T3 to a representative point of each place for i := 1 to #plcs do pt := Eltseq(RepresentativePoint(plcs[i])); assert pt[3] eq 1; // affine point T3pt := SequenceToMultiset(T3(pt[1], pt[2])); F := Universe(Eltseq(pt)); // field of definition // find Frobenius orbits over F orbs := []; while not IsEmpty(T3pt) do t := Rep(T3pt); orb := {Universe(T3pt)| }; t1 := t; repeat Exclude(~T3pt, t1); Include(~orb, t1); t1 := ; until t1 eq t; Append(~orbs, orb); end while; Dpt := Parent(D)!0; for orb in orbs do tpt := Rep(orb); tpt := [tpt[1], tpt[2], 1]; plc := place_from_point(X(Universe(tpt))!tpt); Dpt +:= ExactQuotient(Degree(F)*#orb, Degree(plc))*plc; end for; assert Degree(Dpt) eq 4*Degree(plcs[i]); // sanity check result +:= mults[i]*Dpt; end for; return result; end function; // For a prime p, return a sequence of non-cuspidal places D of degree d, // modulo the action of the diamond operators, such that // t(D) = (T_3 - 3 <3> - 1) ( - 1) (D) is a principal divisor on X_1(p)/F_2. function surviving_places(p, a, d) printf "\nsurviving_places(%o, %o, %o):\n", p, a, d; // get X_1(p) over F_2 X := X1p_models[p]; F := GF(2, d); // get a list of all elliptic curves over F_{2^d} with a point of order p Es := [e[1] : e in curves[d] | IsDivisibleBy(e[2], p)]; printf " there %o %o elliptic curve%o (mod Frob) with a point of order %o\n", #Es eq 1 select "is" else "are", #Es, #Es eq 1 select "" else "s", p; // take a point of order p and turn the curves into points on X // (since the diamond operators permute the pairs (E, P), with P a point // of order p on E over F, transitively, we get a set of representatives // modulo the action of the diamond operators in this way) points := [X(F)| ]; for E in Es do G, mG := AbelianGroup(E); m := #OrderedGenerators(G); assert IsDivisibleBy(Invariants(G)[m], p); pt := mG(ExactQuotient(#G, p)*G.m); Append(~points, X(F)![x,y] where x, y := EP_to_xy(E, pt)); end for; // convert into places places := [Places(X)| place_from_point(pt) : pt in points]; printf " giving %o place%o on X_1(%o)\n", #places, #places eq 1 select "" else "s", p; // now, for every place, compute t(D) and check whether this is principal result := [Universe(places)| ]; for plc in places do D := diamondop_div(a, 1*plc) - plc; D := T3_div(D) - 3*diamondop_div(3, D) - D; printf " testing divisor for being principal...\n"; time flag := IsPrincipal(D); printf " --> %o\n", flag; if flag then Append(~result, plc); end if; end for; printf " %o place%o survive%o.\n", #result, #result eq 1 select "" else "s", #result eq 1 select "s" else ""; return result, places; end function; // Check that no places of degree <= 7 survive for p = 59, 61, 67, 73 // and a suitable a. // We use a = 3 in all cases. This makes the degrees of the positive // and negative parts of the resulting divisors somewhat smaller. printf "\np = 59 (takes about 40 min): "; a := 3; time for d := 1 to 7 do assert IsEmpty(surviving_places(59, a, d)); end for; // ~ 40 min printf "OK!\n"; printf "\np = 61 (takes about 45 min): "; a := 3; // see above, +- the sixth power of a primitive root time for d := 1 to 7 do assert IsEmpty(surviving_places(61, a, d)); end for; // ~ 45 min printf "OK!\n"; printf "\np = 67 (takes about 3 hours): "; a := 3; time for d := 1 to 7 do assert IsEmpty(surviving_places(67, a, d)); end for; // ~ 3 h printf "OK!\n"; printf "\np = 73 (takes about 7 hours): "; a := 3; time for d := 1 to 7 do assert IsEmpty(surviving_places(73, a, d)); end for; // ~ 7 h printf "OK!\n"; //========================================================================= printf "\n=============================================================================\n"; printf "We prove Lemmas 8.8 and 8.9.\n"; printf "=============================================================================\n\n"; // Deal with p = 37 and d = 6, 7. // First show that only one diamond orbit of degree-6 places on X_1(37) // satisfies the principality criterion and no degree-7 place does. // We can take any a here, since the only positive-rank component // of J_1(37) comes from J_0(37) (see below). printf "Show that only one diamond orbit of places of degree 6\n"; printf "satisfies the critrion of Lemma 8.6 and that no other place\n"; printf "of degree d <= 7 does...\n"; for d := 1 to 5 do assert IsEmpty(surviving_places(37, 3, d)); end for; time surv, places := surviving_places(37, 3, 6); // ~ 20 sec assert #surv eq 1; time assert IsEmpty(surviving_places(37, 3, 7)); // ~ 5 sec // Verify that the places do not survive when using the "wrong" operator. for plc in places do D := diamondop_div(3, 1*plc) - plc; D := T3_div(D) - diamondop_div(3, D) - 3*D; time assert not IsPrincipal(D); // ~ 5 sec end for; // Precision bound for q-expansions prec := 200; printf "\nNow show that this gives an essentially unique sporadic point\n"; printf "in X_1(37)^(6)(Q) and no non-cuspidal point in X_1(37)^(7)(Q).\n\n"; printf "\nsetting up modular symbols...\n\n"; // Set up modular symbols for Gamma_1(37). M := ModularSymbols(Gamma1(37), 2, 1); MS := CuspidalSubspace(M); g := Dimension(MS); // 40, the genus of X_1(37) // Get a basis of differentials with integral q-expansions. qbas := qIntegralBasis(MS, prec); // Find quadratic relations between them. printf "constructing canonical model of X_1(37) over F_2...\n\n"; P := PolynomialRing(Integers(), g); mons2 := MonomialsOfDegree(P, 2); monsq := [Evaluate(m, qbas) : m in mons2]; mat := Matrix([[Coefficient(s, j) : j in [2..prec]] : s in monsq]); kermat := KernelMatrix(mat); // takes a few seconds pols := [&+[kermat[i,j]*mons2[j] : j in [1..#mons2]] : i in [1..Nrows(kermat)]]; // Check that the space of quadrics vanishing on the image // has the correct dimension. // Note: binom(g+1, 2) is number of monomials of degree 2, // the dimension of L(2*canonical) is 3(g-1). assert #pols eq Binomial(g+1, 2) - 3*(g-1); // Reduce the model of X_H mod 2. PF2 := PolynomialRing(GF(2), g); pols2 := [PF2!pol : pol in pols]; Pr := ProjectiveSpace(PF2); X := Curve(Pr, pols2); printf "find points over F_{2^6}...\n\n"; // We want points over F_{2^6}, so base-change to that field. X_6 := BaseChange(X, GF(2^6)); // Find all points on X over that field. time pts6 := &join{Points(Scheme(X_6, X_6.2-a*X_6.1)) : a in GF(2^6)} join Points(Scheme(X_6, X_6.1)); // about 5 minutes // Partition into orbits under Frobenius. orbits := frob_orbits(pts6); printf "determine factor of winding quotient...\n\n"; // The only positive-rank component is the first space in the decomposition, // corresponding to the elliptic curve of conductor 37 and rank 1. decomp := NewformDecomposition(MS); assert [Dimension(D) : D in decomp] eq [1,1,2,2,4,6,6,18]; S := &+decomp[[3,5,6,7,8]]; // Verify that the associated abelian variety belongs to the // winding quotient and all isogeneous abelian varieties // have odd order torsion. assert forall{M : M in decomp[2..#decomp] | not(IsZero(WindingElementProjection(M)))}; assert IsOdd(TorsionBound(S, 10)); // Verify that T_{17} has eigenvalue 0 on the positive-rank space // and has nonzero eigenvalues mod 2 on S assert IsZero(HeckeOperator(decomp[1], 17)); assert IsOdd(Determinant(IntegralHeckeOperator(S, 17))); printf "set up projection...\n\n"; // Find the corresponding projection phi: X --> P^{dim A - 1}. // Get Z-basis of differentials in terms of q-expansions. qbasS := qIntegralBasis(S, prec); // Write them as linear combinations on the basis of MS. matMS := Matrix([[Coefficient(s, j) : j in [1..prec-1]] : s in qbas]); qbasSinqbas := Matrix([Solution(matMS, Vector([Coefficient(b, j) : j in [1..prec-1]])) : b in qbasS]); // This gives the linear forms over F_2 realizing the projection phi. linforms := [&+[qbasSinqbas[i,j]*PF2.j : j in [1..Ncols(qbasSinqbas)]] : i in [1..Nrows(qbasSinqbas)]]; // Pick representatives of the honest degree 6 points. reps := [Rep(o) : o in orbits | #o eq 6]; printf "verify criterion for d = 6 ...\n\n"; // In principle, we only need to check the conditions on one of the // orbits that is the image of an orbit of degree-6 points on X_1(37) // that survive the principality criterion, // but it is hard to figure out which one it is. // So we do it for all representative points we found. // Evaluate at each point in reps to check that not all are zero. lfeval := [[Evaluate(l, Eltseq(pt)) : l in linforms] : pt in reps]; assert forall{a : a in lfeval | exists{b : b in a | b ne 0}}; // Check for each point that its coordinates generate a 6-dim'l space over F_2. mattr := [Transpose(Matrix([Eltseq(b) : b in a])) : a in lfeval]; assert forall{m : m in mattr | Rank(m) eq 6}; printf "successful!\n\n"; printf "verify criterion for d = 7 ...\n\n"; // We check that the criterion is also satisfied for the sum of // the image of a rational cusp with a divisor of degree 6. // First find the F_2-rational cusps. time pts1 := Points(Scheme(X, X.1)) join Points(Scheme(X, X.2)) join Points(Scheme(X, X.1+X.2)); // ~ 10 sec assert #pts1 eq 18; // exactly the images of the rational cusps // Check non-vanishing condition at the cusps assert forall{pt : pt in pts1 | exists{l : l in linforms | Evaluate(l, Eltseq(pt)) ne 0}}; // Generate rows for the cusps vectors := [Vector(GF(2), [Evaluate(l, Eltseq(pt)) : l in linforms]) : pt in pts1]; // Check the rank condition. assert forall{m : m in mattr, v in vectors | Rank(Matrix(Append(Rows(m), v))) eq 7}; printf "successful!\n\n"; t := Cputime(t0); min := Floor(t/60); sec := t - 60*min; h := Floor(min/60); min := min - 60*h; printf "Total time: %o hours, %o minutes, %o seconds\n\n", h, min, sec;