Moduli-friendly Eisenstein series over the p-adics and the computation of modular Galois representations

We show how our p-adic method to compute Galois representations occurring in the torsion of Jacobians of algebraic curves can be adapted to modular curves. The main ingredient is the use of “moduli-friendly” Eisenstein series introduced by Makdisi, which allow us to evaluate modular forms at p-adic points of modular curves and dispenses us of the need for equations of modular curves and for q-expansion computations in the construction of models of modular Jacobians. The resulting algorithm compares very favourably to our complex-analytic method.


Introduction
In this article, given an integer k and a subgroup of SL 2 (Z) of finite index, we will denote by M k ( ) (resp.S k ( ), E k ( )) the space of modular forms (resp.cusp forms, Eisenstein series) of weight k and level .
Let f = q + n 2 a n q n ∈ S k 1 (N ) be a newform, let l be a finite prime of the number field Q(a n , n 2), and let ρ f,l : Gal(Q/Q) −→ GL 2 (F l ) be the mod l Galois representation attached to f .Our goal in this article is to compute ρ f,l explicitly, that is to say, to find a squarefree polynomial F (x) ∈ Q[x] and a bijection between the roots of F (x) in Q and the nonzero points of the representation space F 2 l such that the Galois action on these roots matches the representation ρ f,l .
(1.0.1)Indeed, these data allow one to determine efficiently the image by ρ f,l of Frobenius elements, thanks to the Dokchitsers' method [1].Following Couveignes's and Edixhoven's ideas [2], we presented in [3] and in [4] a method to compute ρ f,l explicitly.This method is based on complex-analytic geometry and Makdisi's algorithms [5,6] to compute in Jacobians of curves.It is fairly general, but suffers from some limitations, see Sect.6.2 below.
Later on, in [7], we presented another method, based on an adaptation of Makdisi's algorithms to a p-adic setting, to compute Galois representations occurring in the torsion of Jacobians of any (not necessarily modular) algebraic curve given by an explicit model (a plane equation, for example).This method proceeds by computing torsion points over F p , and then lifting them p-adically.
The goal of this article is to present an adaptation of the p-adic method introduced in [7] to the particular case of modular curves, so as to compute representations such as ρ f,l by a p-adic approximation.This method may be roughly summarised as follows: • Let ∈ N be the prime below l, and let N = N , so that the representation ρ f,l occurs (up to twist) in a piece T of the -torsion of the Jacobian J 1 (N ) of the modular curve X 1 (N ).• Fix a prime p 6N , and an unramified extension Q q of Q p large enough that the points of T are defined over it.Let F q be its residue field.• Find sufficiently many pairs (W, P), where W is the Weierstrass equation of an elliptic curve E W over Q with good reduction at p, and P ∈ E W (Q q ) has exact order N .• These pairs represent Q q -points on X 1 (N ).Evaluate a modular form of level N and weight 2 at these points.This makes sense because the Weierstrass equation determines a local trivialisation of the sheaf of modular forms at the corresponding point of X 1 (N ).Repeat this process until procurement of the values of a basis of M 2 1 (N ) at these points.• These values allow us to faithfully represent subspaces of M 2 1 (N ) , and thus to use Makdisi's algorithms to compute in J 1 (N )(Q q ) and J 1 (N )(F q ).• Pick a random point in J 1 (N )(F q ).Use Makdisi's algorithms to multiply it by the appropriate integer to make it -torsion, and then use the action of Frob p to project it onto T .Repeat this process until a generating set of T is obtained.• Lift this generating set p-adically to some chosen accuracy O(p e ), e ∈ N, and then use Makdisi's algorithms to compute all the points of T to this accuracy.• Construct a rational map α ∈ Q J 1 (N ) , evaluate it at the points of T , and form the polynomial 0 =t∈T x − α(t) .This is an approximation mod p e of the polynomial F (x). • Use this approximation to identify F (x) as an element of Q[x].
The output of this method is thus actually the reduction mod p e of a polynomial F (x) as above along with an indexation of the p-adic roots of F (x) by F 2 l \ {0}, where the precision e ∈ N must be high enough to identify F (x) as an element of Q[x].The correctness of this identification may then be confirmed a posteriori by [8].
This method applies to all pairs (f, l) for which there exists a prime p ∈ N satisfying a technical condition (condition (A.5.3)below).We expect this to be the case as long as the "multiplicity one condition" [9, 3.3] is satisfied, which is a very mild requirement [9,Theorem 3.5].
Compared to a direct application of [7], the advantage of this method is that it suppresses the need for plane equations of modular curves, whilst making an extremely limited use of q-expansions (see Remark 4.2.1).
In exchange, this new approach requires evaluating modular forms at p-adic points of modular curves, which is non-trivial since one cannot use q-expansions for this purpose.We overcome this difficulty by using modular forms introduced in [10], which we manage to evaluate p-adically thanks to their transparent interpretation in terms of the moduli problem parametrised by the modular curve.This is why we are able to compute in modular Jacobians without requiring equations for the corresponding modular curve.
More specifically, the techniques introduced in [10] make it possible to compute in the Jacobian of a modular curve of level N ∈ N by using only the coordinates of the N -torsion points of a single elliptic curve E as input data.Moreover, these techniques are completely algebraic, making them compatible with base-change and thus usable over p-adic fields (where p 6N ), finite fields (of characteristic coprime to 6N ), and intermediate objects such as Z/p e Z with arbitrary p-adic precision e ∈ N.
We thus manage to evaluate modular forms at p-adic points of modular curves, to construct p-adic models of modular Jacobians, and hence to adapt [7] into a method to compute explicitly modular Galois representations ρ f,l since these representations occur in the torsion of modular Jacobians.
This article is organised as follows.We begin by gathering in Sect. 2 arithmetic results about modular curves, and in particular, about their cusps and the Galois action on them, most of which are probably well-known to experts but are unfortunately scattered across the literature.Next, we recall in more detail the ideas behind our p-adic method [7] in Sect.3, and explain how to adapt this method to modular curves, provided that we are able to evaluate modular forms at p-adic points of modular curves.Then, in Sect.4, we recall the definition and some of the properties of Makdisi's moduli-friendly Eisenstein series, and explain how they can be used to evaluate modular forms at p-adic points of modular curves.In order to compute modular Galois representations, it then remains to construct "evaluation maps" from modular Jacobians to A 1 , which we do in Sect. 5. We then demonstrate in Sect.6 that our implementation of the p-adic method presented in this article using [11]'s C library outperforms our [12] implementation of the complexanalytic method by a factor ranging from 20 to 100, meaning that computations of moderately "small" modular Galois representations now take minutes instead of hours, and we explain this difference of performance.Finally, in Sect.6.3, we detail various algorithmic optimisations of the method presented in this article, which are required to obtain good timings but whose presentation we had delayed to this final Section for the sake of clarity.

Reminders on modular curves and their cusps
We begin by reviewing some facts which we will need later about modular curves and modular Galois representations.

Classical congruence subgroups and their moduli problems
Let N ∈ N, and define as usual and more generally, given a subgroup H (Z/N Z) × , Denote the corresponding modular curves by X(N ), X 0 (N ), X 1 (N ), and X H (N ).Note that since −1 ∈ SL 2 (Z) acts trivially on the upper-half plane, we have X H (N ) = X H,−1 (N ) where H, −1 denotes the subgroup of (Z/N Z) × generated by H and −1, so we will restrict our attention to the subgroups H (Z/N Z) × that contain −1.
Let us briefly recall the interest of these modular curves and use this occasion to fix some notation and conventions, which we will use throughout the rest of this article.Informally speaking, the curve X(N ) parametrises the pairs (E, β) up to isomorphism, where E is an elliptic curve and is a group isomorphism mapping the standard basis [1,0], [0, 1] of (Z/N Z) 2 to points P, Q ∈ E[N ] such that the Weil paring e N (P, Q) is a fixed primitive N -th root of 1.
We pause here to mention that we normalise the Weil-pairing as in [13, 7.4]: For any ω 1 , ω 2 ∈ C such that Im ω 1 ω 2 > 0, we have [13, 1.3] e N (ω 1 /N, ω 2 /N ) = e 2π i/N (2.1.1)on the elliptic curve C/(Zω 1 ⊕ Zω 2 ); beware that some authors (and [11]) use the opposite normalisation, namely e N (ω 2 /N, ω 1 /N ) = e 2π i/N .This choice of normalisation will matter later (see Theorem 4.2.2 below).We will always view the elements of (Z/N Z) 2 as row vectors; we then have a left action of SL 2 (Z/N Z) on X(N ) defined by where γ ∈ SL 2 (Z/N Z) and β γ is the isomorphism between (Z/N Z) 2 and E[N ] taking v ∈ (Z/N Z) 2 to β(vγ ).It follows that each fibre of the projection map X(N ) → X(1) at an elliptic curve E having no automorphisms other than ±1 is a torsor under SL 2 (Z/N Z)/±1.
Similarly, the curve X 1 (N ) parametrises isomorphism classes of pairs (E, Q), where E is an elliptic curve and Q ∈ E is a point of exact order N ; more generally, X H (N ) parametrises isomorphism classes of pairs (E, H • Q) of elliptic curves equipped with a point of order N up to multiplication by H.The projection map from X(N ) to X H (N ) is given by so that given γ , γ ∈ SL 2 (Z/N Z), the points (E, β γ ) and (E, β γ ) of X(N ) project to the same point of X H (N ) iff γ and γ have the same bottom row up to scaling by H (remember that we are assuming that −1 ∈ H).It follows that given an elliptic curve E such that Aut(E) = {±1} and an isomorphism β : where γ denotes any element of SL 2 (Z/N Z) whose bottom row is [c, d].
The group 1 (N ) is normal in H (N ), and we have the isomorphism Let μ N ⊂ Q denote the group of N -th roots of 1, and identify the Galois group Gal( The moduli interpretation of X(N ) (resp. of X 1 (N ), X 0 (N ), and more generally, X H (N )) makes sense over Q(μ N ) (resp.over Q), so this curve admits a model over Q(μ N ) (resp.over Q) which is compatible with this moduli interpretation.In particular, the diamond operators are defined over Q.For what follows, we must describe such a model precisely.
As in [14, 6.2], consider the function field where j is the j-invariant and for each nonzero v = (c v , d v ) ∈ (Z/N Z) 2 , the modular function f v 0 is defined on the upper-half plane by where and ℘ τ is the Weierstrass ℘ function attached to the lattice spanned by τ and 1.Then we have [14, 6.2] As ℘ τ (z) is an even function of z for all τ , we have In view of (2.1.2),given a subgroup H (Z/N Z) × containing −1, we may thus set [13, 7.7] hence fixing a model for X H (N ) over Q for each H, and in particular, for X 0 (N ) and X 1 (N ).
In particular, for N 5, X 1 (N ) is the moduli space for elliptic curves E equipped with a torsion point of exact order N , or equivalently with an embedding Z/N Z → E[N ], as opposed to an embedding μ N → E[N ]; see [15, 9.3] and Example 2.3.4 below.

Modular Galois representations in modular Jacobians
Let us define a newform of level H (N ) as a newform of level 1 (N ) whose nebenty- denote the finite set of newforms of weight k and level H (N ).This set is acted on by Gal(Q/Q) via the coefficients of q-expansions at the cusp ∞: if f = ∞ n=1 a n (f )q n has nebentypus ε, then the coefficients {a n (f )} n are algebraic integers, and for τ ∈ Gal(Q/Q) we have τ f ∈ N k ( H (N )) with q-expansion (τ f )(q) = n τ (a n (f ))q n and nebentypus τ • ε.Denote the Galois orbit of f by [f ], and let be the set of such Galois orbits.Then the Jacobian J H (N ) of X H (N ) decomposes up to isogeny over Q as where σ 0 (n) = d|n 1 is the number of divisors of n, H M (Z/MZ) × denotes the image of H in (Z/MZ) × , and for each [f ], the Abelian variety can be thought of as the piece of J H (N ) where the Hecke algebra T of weight 2 and level H (N ) acts with the eigenvalue system of f ; more precisely, A [f ] is defined by where Let f be a newform of weight k, level N , and nebentypus ε f ; let K f be the number field Q(a n (f ) | n 2), which contains the values of ε f , let l be a finite prime of K f , and finally, let ∈ N be the prime below l.Suppose we wish to compute the mod l Galois representation ρ f,l attached to f .By [16, 2.1], this representation is also attached to a form whose level is coprime to , so we assume that N from now on.Similarly, by Theorem 2.7 of [9], up to twist by the mod cyclotomic character, this representation is also attached to a form of the same level and of weight comprised between 2 and + 1, so we suppose 2 k + 1 from now on.
If k = 2, then ρ f,l occurs in the -torsion of the Jacobian J 1 (N ) of X 1 (N ).Else, recall [9, p.178] that there exists an eigenform f 2 of weight 2 but level N and a prime l 2 of K f 2 above such that We thus set so in view of (2.2.1), ρ f,l occurs in J H (N ) provided that H ker ε 2 , where ε 2 is the nebentypus of f 2 .Taking H = ker ε 2 , we thus get a modular curve whose Jacobian contains ρ f,l , but whose genus is (hopefully) smaller than that of X 1 (N ), making explicit computations with it more efficient.This is our reason for introducing the modular curves X H (N ); this idea originates from [17, 4.1].More explicitly, (2.2.2) implies that (2.2.4) Remark 2.2.5 Naturally, in many cases, H is a very small subgroup of (Z/N Z) × , so the genus of X H (N ) is the same as, or not much smaller than, that of X 1 (N ).However, there are also cases when the genus of X H (N ) is dramatically smaller than that of X 1 (N ), which makes it possible to compute Galois representations that would otherwise be out of reach; see [4] for some examples.

Remark 2.2.6
In principle, it would be even better to compute ρ f,l directly in the Abelian variety , but the method that we are adapting only works with Jacobians.

Moduli interpretation and Galois action
Let n ∈ N. Recall that the Néron n-gon is the variety C n obtained by gluing n copies of P 1 indexed by Z/nZ by the relation for all i ∈ Z/nZ.Its regular locus is thus the group variety A morphism between such n-gons is an algebraic variety morphism inducing a group variety morphism on the regular locus.Whereas the non-cuspidal points of the modular curve X 1 (N ) correspond to isomorphism classes of pairs formed by an elliptic curve E and a torsion point of E of exact order N , the cusps of X 1 (N ) correspond to isomorphism classes of pairs formed by a Néron n-gon for some n | N equipped with a torsion point of C reg n of exact order N whose multiples meet every component, see [15, 9.3] or [18,II].Such a pair is thus of the form where ζ ∈ μ N and i ∈ (Z/nZ) × are such that the gcd of N , the order of ζ , and the order of i is 1; in particular, it is defined over Q(ζ ) ⊆ Q(μ N ).In order to understand the cusps of X 1 (N ), and how they are permuted by Gal(Q/Q), we must therefore classify such pairs up to isomorphism, and in particular, determine the automorphisms of C n .
First of all, we have canonically Therefore, Finally, an automorphism of C reg n extends to an automorphism of C n iff it respects the gluing condition (∞, i) ∼ (0, i + 1), which translates into the condition j = m.Therefore, and elementary arithmetic considerations, which we omit here for brevity, then show that we have a bijection where for brevity we have written (c, N ) for gcd(c, N ), and where ζ N is a fixed primitive Nth root of unity.We may thus represent the cusps of X 1 (N ) by such pairs (c, d) up to negation.
The advantage of this representation is that the Galois action on the cusps is then transparent.Indeed, (2.3.2) confirms that the cusps are all defined over Q(μ N ), and shows that for each x ∈ (Z/N Z) × , we have where Example 2.3.4 The cusp ∞ = 1/0 is represented by (c = 0, d = 1).In particular, (2.3.3)shows that this cusp is fixed by σ x only if x = ±1, so its field of definition is ).This can be visualised by noticing that it corresponds under (2.3.2) to the pair C 1 , (ζ N , 0) , and by applying (2.3.1) to n = N /(c, N ) = 1, which shows that this pair is isomorphic (by negation) to C 1 , (ζ −1 N , 0) but not to any other of its Galois conjugates.
On the contrary, the cusp 0 = 0/1 is represented by (c = 1, d = 0), and is thus defined over Q; indeed, it corresponds to the pair C N , (1, 1) , which is clearly defined over Q.

Moduli interpretation of the cusp 0
As explained in the previous Subsection, in this article, we will actually not work with X 1 (N ), but rather with X H (N ) where H is a subgroup of (Z/N Z) × containing −1.Fortunately, translating the above results to the case of the cusps of X H (N ) presents no difficulty.Indeed, going down from X 1 (N ) to X H (N ) amounts to identifying

Widths
In order to consider q-expansions at various cusps, we must determine the width of these cusps.
Let s be a cusp of X H (N ), and let again M s ∈ SL 2 (Z) be such that M s • ∞ = s.The width of s is then by definition the smallest positive integer w such that

Rationality of q-expansions
In order to compute modular Galois representations, we will need to construct rational maps J H (N ) A 1 .As we will see in Sect.5, one way to do so involves looking at the qexpansion coefficients a n (f ) of some modular forms f at some cusp; but it is of course fundamental for our purpose that the dependency of the a n (f ) on f be Galois-equivariant.This is unfortunately not the case at every cusp; for instance, it is not the case of the cusp ∞ of X 1 (N ) for N > 6 since this cusp is not even defined over Q.
In fact, thinking about the rationality of the q-expansion in terms of the cusp alone is wrong.Indeed, given a function f ∈ Q X H (N ) (or more generally, a modular form) and a cusp s of X H (N ) of width w, it is tempting to define "the" q-expansion of f at s as the expansion of f M s at ∞ in terms of q w = e 2π iτ /w , where M s ∈ SL 2 (Z) satisfies M s • ∞ = s.However, this definition does not make actual sense if w > 1.Indeed, the matrices M ∈ SL 2 (Z) satisfying M • ∞ = s are precisely those of the form M = ±M s 1 x 0 1 where x ∈ Z, and while the ± sign does not matter, different values of x yield different q w -expansions of f M at ∞; more precisely, we have for all n ∈ Z.However, this still shows that the coefficient a 0 (f M s ) and the order of vanishing of f do not depend on M s , but only on s.
We are thus led to the following definition: Definition 2.3.14 We say that a matrix M ∈ SL 2 (Z/N Z) yields rational q-expansions if the map This condition is equivalent to the requirement that the q-expansion of f M have rational coefficients whenever f is defined over Q.
Remark 2.3.15Technically, these expansions are q w -expansions, where w ∈ N is the width of the cusp M • ∞.One way to circumvent this technicality would be to talk about q N -expansions, since w | N always.However, this does not impact our discussion about the Galois-equivariance of the coefficients, so for convenience, we will persist in this abuse of language in the rest of this Section.
We must therefore determine which M ∈ SL 2 (Z/N Z) yield rational q-expansions.Recall from [14, 6.2] that every element f ∈ F N has a (possibly Laurent) q N -expansion , and that we have the relation for all x ∈ (Z/N Z) × and n ∈ Z, where σ x ∈ Gal(Q(μ N )/Q) is as in (2.1.6).From (2.1.8),we deduce that for all M ∈ SL 2 (Z/N Z), This criterion allows us to determine explicitly for which cusps s of X H (N ) there exists M ∈ SL 2 (Z/N Z) such that M • ∞ = s and that M yields rational q-expansions, and to find such an M if one exists.

Remark 2.3.18
There is always at least one such cusp; namely, we can take s = 0 and M = 0 −1 3 Constructing p-adic models of modular Jacobians

Computing p-adically in Jacobians
Let us summarise the kind of data that we need to describe a curve in whose Jacobian we want to compute Galois representations p-adically by using our methods presented in [7], which are themselves based on Makdisi's algorithms [5,6].After this, we will expose in Sect.3.2 below how to obtain these data in the specific case of the modular curves X H (N ) which we use to compute modular Galois representations.
Definition 3.1.1Let C be a projective, geometrically non-singular curve of genus g defined over Q.Let p ∈ N be a prime number at which C has good reduction, a ∈ N an integer, q = p a , Q q the unramified extension of Q p of degree a, Z q its ring of integers, F q its residue field, and finally let e ∈ N. A p-adic Makdisi model of C with residue degree a and p-adic accuracy O(p e ) consists of: • A choice of a line bundle L on C whose degree d 0 = deg L satisfies which reduce to pairwise distinct points of C(F q ), and which are globally invariant under Frob p , • The permutation induced by Frob p on the points generally, over Q p ) at each of the points P i , so that we have a Galois-equivariant concept of value of a global section of L at each P i ; , where d = dim H 0 (C, L) = d 0 + 1 − g and the v j form some Q q -basis of H 0 (C, L) such that the values v j (P i ) under the t i lie in Z q and which remains an that is to say, the numerator of the Zeta function of C/F p reversed so that it is of the form Similarly, the bound (3.1.2) ensures that we avoid complications stemming from dealing with Riemann-Roch spaces attached to divisors of low degree so that Makdisi's algorithms to compute in Jacobians are valid.For instance, it ensures that the multiplication map is surjective, so we can compute the global sections of L ⊗n , n 5 from the datum of the matrix V .
Remark 3.1.6Note that in particular, such a p-adic Makdisi model of C does not include an explicit model for C.This is because (3.1.2) ensures that L is very ample and thus defines a projective embedding of C, whose equations could be read off the kernel of multiplication maps such as (3.1.5).But of course, to construct such a p-adic Makdisi model of C, some explicit data about C must be known so as to write down the matrix V .In this article, C will be a modular curve, and the v j will be modular forms constructed in Sect.4.4.
By Riemann-Roch and assumption (3.1.2),every point x ∈ J = Pic 0 (C) is represented by the line bundle L(−D x ) for some (non-unique) effective divisor D x on C of degree d 0 = deg L. A point x ∈ J (Q q ) may thus be represented by the matrix where the w j form a Q q -basis of the space of global sections of L ⊗2 (−D x ) chosen so that the values w j (P i ) under the local trivialisations t i lie in Z q and that the w j still determine a basis of this section space over F q .
As explained in [7], this mode of representation of points of J allow us to perform padic computations in J with a p-adic Makdisi model of C thanks to Makdisi's algorithms without loss of p-adic accuracy, which leads to a method to compute the representations of Gal(Q/Q) afforded in the torsion of the Jacobian J of C.
More precisely, let p be a prime, let T ⊆ J [ ] be a Galois-submodule, let be the mod Galois representation afforded by T , and denote by the characteristic polynomial of the Frobenius Frob p acting on T .If a = [F q : F p ] has been chosen so that Gal(Q q /Q q ) acts trivially on T, (3.1.8)so that the points of T are defined over Q q , and if χ T (x) is coprime mod with its cofactor L p (x)/χ T (x), (3.1.9)so that the datum of χ T (x) determines the submodule T ⊆ J [ ] non-ambiguously, then we can compute ρ T as follows: 1. Determine the order N = #J(F q ) as N = Res(L p (x), x a − 1), and factor it as 2. Take a random point x ∈ J(F q ), multiply it by M to get an -power torsion point, then repeatedly by to get an -torsion point, and finally apply Lp χ (Frob p ) to project it onto T , 3. Repeat this process to obtain an F -basis of T over F q , 4. Lift this basis of T to J(Z q /p e )[ ], 5. Use Makdisi's algorithms over Z q /p e to compute all the points of T over Z q /p e , by forming all linear combinations from this basis, 6. Construct a rational map α : J A 1 defined over Q, , and identify it as an element of Q[x].

Strategy 3.1.10 p-adic computation of Galois representations found in Jacobians
Indeed, if α is sufficiently generic to be injective on T , then the values α(t) are permuted by Gal(Q/Q) in a way matching the points of T , so the polynomial F (x) satisfies (1.0.1).

Remark 3.1.11
We actually get a mod p e approximation of F (x), so we need the accuracy parameter e to be large enough to identify F (x) ∈ Q[x].Naturally, higher values of e will also work, but slow the computation down.As in [7], we do not have a clear recipe for the ideal value of e, and we proceed mostly by trial-and-error.Therefore, the results which we get are not rigorously proved to be correct; however, in the case of Galois representations attached to modular forms, they can be rigorously certified using the methods presented in [8].In what follows, we will not concern ourselves with this aspect anymore and assume that the value of e has been set somehow.
Remark 3.1.12 The reason why we record how the Frobenius Frob p permutes the points P 1 , • • • , P n Z in the p-adic Makdisi model is that this allows us to apply Frob p to points of J , see Sect.A.1 below for details.

Overview of the construction of p-adic Makdisi models of modular Jacobians
Suppose we want to compute the mod l representation ρ f,l attached to a newform f of weight k, level N , and nebentypus ε f by Strategy 3.1.10.As explained in Sect.2.2, this representation is found up to twist in the -torsion of the Jacobian J H (N ) of X H (N ), where N ∈ N is defined by (2.2.3) and H (Z/N Z) × is defined by (2.2.4), so we must construct a p-adic Makdisi model of X H (N ).
We are going to assume for now that the value of the prime p has been fixed; we will then explain in Sect.3.3 below how to optimise this value.
For simplicity, we make the following assumptions: • X H (N ) has at least 3 cusps, • p does not divide 6, nor , nor N , nor the order of H (Z/N Z) × .
We require p N since p-adic Makdisi models require p to be a prime of good reduction, p = since our method [7] relies on the -torsion being étale at p, p 6 to ensure the validity of Makdisi's construction of the moduli-friendly modular forms which we detail in Sect.4, and p #H so as to facilitate our use of Makdisi's construction (see Sect. 4.4).On the contrary, the assumption on the number on cusps is not an essential one (see Remark 3.2.2below).
The value of p imposes the optimal residue degree a = [Q q : Q p ] to use to compute ρ f,l as the order of ρ f,l (Frob p ).This order can usually be determined explicitly, and in any case bounded, from the datum of the characteristic polynomial of ρ f,l (Frob p ), as explained in [7,Proposition 6.1].This requires evaluating mod l the quantities a p (f ) and ε f (p), which we do by classical methods.
In order to construct a p-adic Makdisi model for X H (N ), we need, in particular, to determine the local L factor of X H (N ) at p.For this, we suppose that we can compute the set of Galois orbits of mod N Dirichlet characters where n (t) ∈ Z[t] denotes the n-th cyclotomic polynomial, and that for each such orbit, we can compute the matrix of the Hecke operator T p with respect to some Q[t]/ ord χ (t)basis of the space of cusp forms S 2 (N , χ); for instance, this is possible using [11].Then, in view of the decomposition (2.2.1), we have In particular, we thus recover the genus g of X H (N ) as g = 1 2 deg L p (x); although this is certainly not the most direct way to determine g, this dispenses us from the need to compute it separately.
We also need to pick a line bundle L on X H (N ).It is natural to choose a line bundle whose sections are modular forms of level H (N ); we choose L so that its sections are the modular forms (not just cusp forms) of weight 2. The degree of L is 2g − 2 + ν ∞ , where ν ∞ is the number of cusps of X H (N ); the assumption ν ∞ 3 that we have made above thus ensures that the requirement (3.1.2) is met.Remark 3.2.2Even if our modular curve has fewer than 3 cusps, we can still apply the same construction, by choosing L so that its sections are modular forms of some higher weight, thus ensuring that deg L is large enough.However, in almost all the cases that will be relevant to us in this article, ν ∞ is much larger than 3, so we assumed that ν ∞ 3 for the simplicity of the exposition.
Next, we need to fix sufficiently many such points P i ∈ X H (N ) at which to evaluate (under some local trivialisations) a basis of global sections of L. Since we wish to avoid the use of a plane model for X H (N ), we will rely on the moduli interpretation of X H (N ) detailed in Sect.2.1, thus viewing points as pairs (E, H • Q) where E is an elliptic curve and . Furthermore, instead of using a different elliptic curve E for each point, we follow Makdisi's brilliant idea [10,Theorem 5.4] to fix the curve E, and to consider the points on X H (N ) corresponding to the various Q ∈ E[N ]; in other words, we will work in a fibre of the projection map π : X H (N ) −→ X (1).
We pick an elliptic curve E defined over Q, having good reduction at p, whose j-invariant is neither 0 nor 1728 mod p so as to avoid the ramification locus of π, and whose N -torsion is defined over Q q .We will show in Sect.3.3 that such curves E always exist provided that p is not too small, and we will explain in detail how to find efficiently such a curve E in Sect.A.4.
As j(E) ≡ 0, 1728 mod p, the number of points in the fibre of π above E is then equal to the degree d of π.By [13, 3.1.1],the genus of X H (N ) is where ν 2 (resp.ν 3 ) denotes the number of elliptic points of X H (N ) of order 2 (resp.3).It follows that which shows that the lower bound (3.1.3)on the number of points of X H (N ) at which we evaluate the sections of L is satisfied; we are thus in good shape to construct a valid Makdisi p-adic model of X H (N ).We also fix arbitrarily an isomorphism β : . By (2.1.3),the points at which we evaluate our forms, that is to say, the points on the fibre of the projection X H (N ) → X(1) at E, may then be identified with the primitive vectors of (Z/N Z) 2 up to scaling by H.In particular, if ∈ GL 2 (Z/N Z) is the matrix describing the action of the Frobenius Frob p on E[N ] with respect to β, then the image by Frob p of the point of the fibre corresponding to the primitive vector v ∈ (Z/N Z) 2 is the point of the fibre corresponding to the vector v •( t ).We can therefore determine the permutation induced by Frob p on the fibre, provided that we have computed the matrix .We explain in detail how we determine in Sect.A.4 below.
To complete the construction of our Makdisi model, the only remaining step is to write down the matrix V .For this, we need to evaluate under some local trivialisations a basis of global sections of L, which are modular forms, at the points of X H (N ) in the fibre of π corresponding to E. This requires ingredients that will be introduced in Sect.4; we will actually show in Sect.4.4 how these values can be read off the coordinates in Q q of the N -torsion points of E. Remark 3.2.3In order to compute Galois representations by Strategy 3.1.10,we need to generate torsion points in the Jacobian over F q , and for this, we must generate random points over F q .In Makdisi's algorithms, a point on the Jacobian is represented by a subspace of a fixed Riemann-Roch space defined by vanishing conditions at an effective divisor on the curve.Bruin [19,Algorithm 3.7] presents a sophisticated method to generate uniformly distributed random points on the Jacobian in this framework.Still, as explained in [7, 6.2.1], we use a much cruder (and faster) approach, which in the case of modular curves amounts to constructing by linear algebra subspaces of modular form spaces consisting of forms that vanish at some of the points P i ∈ X H (N ) at which we evaluate our forms.Since these points are in a rather special configuration, namely as they all lie on the same fibre of the projection to X(1), it may happen that the random points of J H (N )(F q ) obtained this way are so poorly distributed than they generate a subgroup with so little -torsion that it does not allow us to generate the representation space, making the computation of the Galois representation stall at stage 2 of Strategy 3.1.10.Fortunately, this seems rare in practice; in fact, in most of the cases that we have encountered, switching to another elliptic curve E suffices to solve this issue.Another workaround consists in using not one but several elliptic curves E, so as to allow ourselves to work with divisors supported on several fibres of the projection to X(1); and if this also fails, then we can fall back to Bruin's method.

The choice of the prime p
We now explain how to optimise the value of the prime p.
As explained in the previous Subsection, the degree a = [Q q : Q p ] of our Makdisi model must satisfy two constraints, which both depend on p: first, there must exist an elliptic curve as above having its N -torsion defined over Q q , which imposes by the Weil pairing and by the Hasse bound; and second, a must be a multiple of the order of ρ f,l (Frob p ). Naturally, the smaller a, the more efficient the computations will be (bearing in mind the remarks made in [7, 6.4]), so it is advantageous to try many values of p, and to select the one resulting in a being as small as possible.Conditions (3.3.1) and (3.3.2) are necessary, but not sufficient, for there to exist a suitable elliptic curve E, or in other words, for X(N ) to admit an F q -rational point away from the locus j ∈ {∞, 0, 1728}; for instance, there exists no elliptic curve over F 13 having j / ∈ {0, 1728} and full 4-torsion over F 13 .Determining necessary and sufficient conditions on p and a in terms of N is an interesting problem, that could probably be solved by examining the Zeta function of the modular curve X(N ).We have chosen not to go this way; instead, we simply try random Weierstrass equations until we find a curve having full N -torsion over F q and j / ∈ {0, 1728}, and give up and increase p if no such curve is found after a certain number of attempts.Indeed, Cebotarev ensures the existence of infinitely many values of p for which ρ f,l (Frob p ) has the same small order, and the Weil bounds applied to X(N ) show that elliptic curves having full N -torsion over Q q will exist for p large enough.The value of p must be such that χ p (x) is coprime mod l with its cofactor in the Lfactor L p (x) of X H (N ) at p, so that we can isolate the subspace of J H (N )[ ] affording ρ f,l thanks to the action of Frob p , see (3.1.9).This extra constraint does not usually preclude a lot of values of p.However, even in the rare occurrence when no such p can be found, we can still compute ρ f,l by using a variant of Strategy 3.1.10which bypasses (3.1.9)and whose presentation we defer to Sect.A.5.
A reasonable strategy is thus to determine in parallel χ p (x) and L p (x) for all p not dividing 6 N #H up to some bound B, and to retain a value of p satisfying the above constraints and leading to a degree a which is as small as possible.The value of B depends on how fast we can determine χ p (x), which requires evaluating a p (f ) mod l, and L p (x), which by (3.2.1) involves computing the action of the Hecke operator T p on S 2 H (N ) , and how insistent we are to minimise the degree a; in practice, we use B = 100 or 1000, see the examples in Sect.6 below.

Makdisi's construction
Makdisi introduces in [10] the concept of moduli-friendly modular forms, which are forms whose value at a point of the modular curve can be easily read off the representation of this point as an elliptic curve equipped with some appropriate level structure.We are now going to review this construction, after what we will explain in detail in Sect.4.4 how to use it to write down the matrix V (see Definition 3.1.1)so as to complete our construction of p-adic Makdisi models of modular curves X H (N ) without resorting to explicit plane models, since we had left this question open in the previous Section.
From now on, we think of modular forms "à la Katz" [20]; in other words, we view a modular form f of weight k and level (N ) over a ring R in which N is invertible as a function on the set of isomorphism classes of triples where E is an elliptic curve over R, ω is a generator of the sheaf of regular relative differentials on E/R, and β is what [20, 2.0.3]calls a naïve level N structure on E, that is to say, an isomorphism of group schemes over R, and satisfying the homogeneity condition for all λ ∈ R × as well as some extra compatibility conditions (namely commutation with base change, see [20, 2.1] for details).
In fact, we will restrict ourselves to elliptic curves defined by short Weierstrass equations (W) : By assigning to such an equation the differential ω W = dx/2y, it then makes sense to talk about the value f (W, β) of f at the pair (W, β); in other words, choosing a short Weierstrass model for E yields a local trivialisation of the sheaf of modular forms at the corresponding point.Fix a level N 3, and let R be a ring in which 6N is invertible.In [10], Makdisi constructs Eisenstein forms f v 1 of weight k = 1 and level (N ) over R indexed by non-zero vectors v ∈ (Z/N Z) 2 and which enjoy the following particularly nice properties: where W is a short Weierstrass equation defining an elliptic curve E over R and β is a naïve level N structure on E, the value agrees with the slope of the line joining the aligned points β(v 1 ), β(v 2 ), and β(v 3 ) on the model of E defined by W (to be interpreted as the slope of a tangent line in the case where some of in other words, as [10] puts it, it "misses" precisely the cusp forms of weight 1.
When P, Q ∈ E are points such that neither P, nor Q, nor P + Q are at infinity, ( denote by λ P,Q the slope of the line joining P and Q on the model W of E if P = Q, and the slope of the tangent line of E at P if P = Q; this slope is well-defined as (4.1.3)ensures that this line is not vertical.The first property can be summarised by even in the case where v 1 , v 2 and −v 1 − v 2 are not distinct.Makdisi shows that this relation can be inverted so as to read the value of the f v 1 off the slopes of the lines joining the N -torsion points of E, for instance as where w is any vector of (Z/N Z) 2 whose span intersects trivially that of v so that the slope λ β(v),β(v+xw) is well-defined for all x (however, see Proposition 4.3.1 below for a more efficient method).Thus the f v 1 are indeed moduli-friendly modular forms.By combining this observation with the second property, we thus get the neat statement that apart from cusp forms of weight 1, any modular form can (in principle) be expressed as a polynomial in the moduli-friendly forms f v 1 , and therefore evaluated at a pair (W, β).

The q-expansion of f v 1
We briefly pause our construction of Makdisi models of modular curves to determine the q-expansion of the forms f v 1 .
Remark 4.2.1Although we strive not to rely on q-expansions, we will actually make a limited use of them in order to optimise our Makdisi model in Sect.6.5.We could still construct non-optimised Makdisi models of modular curves while completely dispensing of q-expansions; however, to compute modular Galois representations by the methods presented in this article, we must also construct rational maps from J H (N ) to A 1 , and this construction requires a limited use of q-expansions, see Sect. 5 below.
By [10, corollary 3.13], f v 1 (W, β) is proportional to the Weierstrass Zeta function of the elliptic curve defined by (W) evaluated at the N -torsion point β(v), and therefore to the modular form denoted by g v 1 in [13, 4.8].Combining the formulas found in [13], especially in Sect.4.8 thereof, we then find after some computations whose details we omit the following formula. where and for all n 1, Remark 4.2.3It should also be possible to derive these formulas by evaluating the modulifriendly form f v 1 on the Tate curve.Unfortunately, although this would be more in the spirit of this article, this seems to lead to more difficult computations.Remark 4.2.4These formulas allow us to determine the coefficients a n for n up to some bound B in quasi-linear time in B. Since Theorem 4.1.2guarantees that apart from cusp forms of weight 1, every modular form over a congruence subgroup of level N is expressible as a polynomial in the f v 1 , we can thus compute the first B coefficients of the q-expansion of any such form in quasi-linear time in B thanks to fast series arithmetic.This is quasioptimal, and faster than other methods such as modular symbols whose complexity is at least quadratic in B; furthermore, by nature, this approach is well-suited to the Chinese remainder strategy, which involves the computation of the desired result modulo several primes and is the key to many fast algorithms in computer algebra [21,5].However, the bad complexity of this approach with respect to the level N makes it uncompetitive with standard methods unless B is very large.

Efficient evaluation of f v 1
In order to construct p-adic Makdisi models of modular Jacobians, we will use the modular forms f v 1 over rings of the form R = Z q /p e where p 6N , in which the most computationally expensive of the four basic operations is by far division.Using (4.1.5)to evaluate such a form f v 1 at a pair (W, β) requires us to determine N slopes of lines joining points of E[N ], or of tangent lines at such points.Evaluating such a slope requires one division in R, even in the case of a tangent line, since differentiating (4.1.1)shows that the slope of the tangent to E at the point P = (x, y) is λ P,P = 3x 2 +A 2y .In total, each evaluation of f v 1 by (4.1.5)therefore requires N divisions in R, which may be prohibitively costly if e and N are large.We will therefore use a better approach, which was hinted at by Makdisi in [10, 3.14] and allows one to evaluate f v 1 by evaluating only O(log N ) slopes.This approach, which we now describe in detail basing ourselves on explanations provided by Makdiski to the author, is essentially a "double and add" algorithm.Proposition 4.3.1 Let (W, β) be a pair as above, let v ∈ (Z/N Z) 2 be a non-zero vector, and let n | N be the exact order of v, so that P where for brevity we have written P m for 3) is satisfied in every case of (4.3.2),(ii) For all 1 m < n, we have f mv Proof (i) In the case where 1 < m < n is even, neither P m/2 nor P m/2 + P m/2 are at infinity since m/2 < m < n.In the case where 1 < m < n is odd, neither P nor P m−1 nor this is apparent on (4.1.5),and actually follows directly from Makdisi's construction.The formula f mv 1 (W, β) = mf v 1 (W, β) − c m then follows by induction on m < n.Indeed, it obviously holds for m = 1.Suppose now that it holds for all m < m.In the case where m is even, (4.1.4)applied with Similarly, in the case when m is odd, (4.1.4)applied with (iii) Taking m = n − 1 in (ii) and using again the fact that f v 1 is an odd function of v, we obtain Remark 4.3.3 The formula (4.1.5)and the algorithm outlined in Proposition 4.3.1 both demonstrate that the forms l v,w 1 : (W, β) → λ β(v),β(w) span the same R-algebra of modular forms as the f v 1 .However, although the generators l v,w 1 may seem more appealing since they are easier to evaluate than the f v 1 , we shall demonstrate in Remark (4.4.3) below that using the f v 1 results in a better complexity in the construction of p-adic Makdisi models of modular curves, whence our focus on the f v 1 in this Section.

Application to p-adic Makdisi models
In Sect.
where A, B ∈ Z are such that W (not just E) has good reduction at p. Since p N by assumption, Néron-Ogg-Shafarevic [22, Theorem 1] ensures that the coordinates of the N -torsion points of E generate an unramified extension Q q of Q p .Furthermore, our assumption that W has good reduction at p ensures that the coordinates of the nonzero points of E[N ] actually lie in the ring of integers Z q of Q q .
By the second part of Theorem 4.1.2,the space of modular forms of weight 2 and level (N ) is spanned by the products f v 1 f w 1 , where v and w range over the set of nonzero vectors of (Z/N Z) 2 .We use these products f v 1 f w 1 to construct weight-2 forms of level H (N ) by taking traces; namely, we set where These forms do span M 2 ( H (N )), as # H (N ) = N #H is a p-adic unit thanks to our assumption that p #H.We thus proceed as follows to complete the construction of our p-adic Makdisi model of X H (N ).
Having fixed the desired p-adic accuracy O(p e ) for our Makdisi model of X H (N ), we first determine the coordinates in Z q /p e of the points of E[N ] in the model W (we explain how to do this efficiently in Sect.A.4 below).
Recall that we fixed β : (Z/N Z) 2 E[N ] in Sect.3.2.As we explained, the Weierstrass model of E provides us with a normalisation of the differential on E and thus with a local trivialisation of L at the corresponding point, so the value of a modular form at this point is a well-defined quantity.Explicitly, (2.1.3)shows that for all primitive u ∈ (Z/N Z) 2 , the value of f v,w 2,H at the point of X H (N ) represented by W, where U is any element of SL 2 (Z/N Z) whose bottom row is u.We can thus compute this value mod p e thanks to Proposition 4.3.1.
In order to obtain a basis of M 2 H (N ) , which has known dimension we simply successively pick random pairs (v, w) of nonzero vectors of (Z/N Z) 2 , and form for each such pair the vector of values of the form f v,w 2,H at all the points of the fibre, that is to say, the vector of the f v,w 2,H W, H • β(u) where u ranges over the primitive vectors of (Z/N Z) 2 mod H, until the reduction mod p of these vectors has rank d 2 .We then extract a basis, and thus obtain the matrix V for our p-adic Makdisi model.Remark 4.4.2By sticking to the moduli interpretation of modular curves, we have thus managed to obtain a p-adic Makdisi model for X H (N ) without requiring plane equations nor writing down a single q-expansion, merely by looking at the N -torsion of just one elliptic curve E over Q. Besides, this method is straightforward to generalise to modular curves corresponding to any congruence subgroup.It could even be generalised to Shimura curves if an analogue were known for Makdisi's moduli-friendly forms in this context, but sadly this does not seem to be the case at the time of writing.
Remark 4.4.3As mentioned in Remark 4.3.3,we may be tempted to construct our forms of weight 2 by taking products of two forms of the form l v,w 1 : (W, β) → λ β(v),β(w) rather than the f v 1 .The matrix V has size O(g) × O(g), so this would require us to evaluate O(g 2 )O(# H (N )) such forms l v,w 1 , and therefore to perform that many divisions in Z q /p e .In the case where H (N ) = 0 (N ) with N prime, that is O(N 4 ) divisions in R; whereas in the case where H (N ) = 1 (N ) with N prime, that is O(N 5 ) divisions.In contrast, by working with the f v 1 , we can precompute the N × N matrix containing the f v 1 (W, β) for all v ∈ (Z/N Z) 2 , which only requires O(N 2 log N ) divisions with the algorithm outlined in Proposition 4.3.1, after what no further divisions are required to fill in the matrix V .Such a precomputation would not be of any help with the l v,w 1 , since there are O(N 4 ) sets of the form {v, w, −v − w} with v and w in (Z/N Z) 2 , and therefore that many forms l v,w 1 .

Construction of evaluation maps
In order to compute Galois representations following Strategy 3.1.10,the final ingredient that we need is the construction of one or more rational maps α : J H (N ) A 1 defined over Q.
We begin by constructing a rational map where V 2 denotes the space of global sections of L ⊗2 , by the method exposed in [7, 2.2.3].As explained in the same reference, this requires picking two linearly-inequivalent effective divisors E 1 E 2 on X H (N ) of degree d 0 − g, such that we can compute the subspace of V 2 formed by sections that vanish at E i for each i ∈ {1, 2}; besides, these divisors must be defined over Q to ensure that α P is defined over Q.
We choose E 1 and E 2 so that they are supported by cusps, possibly with multiplicities.It is then easy to arrange that E 1 and E 2 are defined over Q since the modular curve X H (N ) tends to have plenty of rational cusps as mentioned in Remark 2.3.9;besides, as explained in Remark 3.1.4,V 2 is spanned by products to two forms f v i ,w i 2,H , whose q-expansions can be determined by (4.4.1) and Theorem 4.2.2, so we can determine the subspaces of V 2 corresponding to these divisors by linear algebra, even if these divisors have multiplicities.
In order to get an A 1 -valued Galois-equivariant map, it remains to construct a rational map A 1 defined over Q.For this, we offer two strategies.

Using q-expansions
The first strategy involves q-expansions; namely, similarly to [3, 3.6], we construct this map as where n 1 , n 2 are nonnegative integers and M 1 , M 2 ∈ SL 2 (Z/N Z) yield rational qexpansions in the sense of Definition 2.3.14.
As noted in Remark 2.3.18,there always exists at least one matrix which yields rational qexpansions, namely M = 0 1 −1 0 .Therefore, this construction applies to every level, for example by taking M 1 = M 2 = M and n 1 = n 2 (lest we get a constant map).In fact, there are always infinitely many possible choices for the parameters n 1 and n 2 , and usually many choices for M 1 and M 2 as well.For instance, we can enumerate the cusps of X H (N ) for which there exists a matrix which yield rational q-expansions, and try all the pairs (M 1 , M 2 ) of (not necessarily distinct) matrices in this list and all integers n 1 , n 2 up to some bound B, typically in the vicinity of B = 5.
We thus obtain several evaluation maps α, some of which may not be injective on the Fsubspace of J H (N )[ ] which affords our Galois representation, and therefore not useful for our purpose; however, in practice, if the bound B is large enough, we always get many injective versions of α, and thus many versions of the polynomial F (x) which describes the Galois representation (see Strategy 3.1.10).We then simply keep the "prettiest" polynomial, for instance that having the smallest arithmetic height.
Remark 5.1.2Recording the q-expansions of forms f i ∈ V 2 forming a basis of V 2 during the creation of the p-adic Makdisi model of X H (N ) makes evaluating (5.1.1)at f ∈ V 2 easy, since we merely have to identify by linear algebra f as a linear combination of the f i from their values at some points of the fibre of X H (N ) → X(1).

Using forms defined over Q
Another strategy consists in fixing a basis (f i ) i of the space V 2 , and in considering the map where i 1 , i 2 dim V 2 are fixed integers.This is in fact the adaptation of the approach that we used in [7, 2.2.3], which applies to any curve, not just modular curves.Again there many choices of pairs (i 1 , i 2 ), so we get many versions of α and of F (x).
However, the basis (f i ) i of V 2 must be made up of forms which are defined over Q for the resulting map to be Galois-equivariant.For this, it is enough to construct a basis of the section space of L formed of forms which are defined over Q, since we can then generate V 2 by taking products of two such forms (see Remark 3.1.4).The sections f v,w 2,H of L introduced in (4.4.1) are, in general, only defined over the cyclotomic field Q(μ N ).However, given a form f defined over Q(μ N ), (2.3.16)applied to the quotient f /f 0 , where f 0 is a form defined over Q and of the same weight as f , combined with the fact that M = 0 1 −1 0 yields rational coefficients, shows that Therefore, for all nonzero v ∈ (Z/N Z) 2 and w ∈ (Z/N Z) 2 and for all y ∈ Z/N Z, the section of L is defined over Q .Furthermore, under the assumption that p φ(N ) (a mild strengthening of the assumption p #H made earlier in Sect.3.2), these sections generate the section space of L over Z p /p e for any e ∈ N, since φ(N )# H (N ) is then a p-adic unit.
Unfortunately, a bit of experimenting with our implementation has revealed that this approach results in polynomials F (x) whose arithmetic height is tremendously larger than those obtained with the approach based on q-expansions, and which therefore require ridiculously high p-adic accuracy to be identifiable as elements of Q[x], see Remark 6.2.1 below.This expresses the fact that products of two forms of the form (5.2.1) do not form "nice" Q-bases of M 4 H (N ) , and the linear algebra used in Sect.6.5 below to optimise the Makdisi model make things even worse.For this reason, in practice we only use the construction of maps α based on q-expansions.

Comparison with the complex-analytic method and results
We conclude be giving some examples to demonstrate the performance of our implementation of the method presented in this article.This implementation benefits from several optimisations, whose presentation we have postponed to Appendix A below for the sake of the exposition.
When analysing these examples, the reader should bear in my mind that the difficulty of the computation of a mod representation is governed by two essential parameters: the genus of the modular curve used in the computation of course, but also the number itself, since the computation must process #(F 2 \{0}) = 2 −1 torsion points, and outputs a polynomial F (x) ∈ Q[x] of degree 2 − 1 whose arithmetic height is likely to grow with , thus requiring more p-adic accuracy.

A form of weight 2 and level 16
The form of [23] label 16.2.e.a is up to Galois-conjugacy the only newform of weight k = 2 and level 1 (16).Since its coefficient field K f = Q(i) is an extension of Q of degree 2, the modular curve X 1 ( 16) is of genus g = 2. Since f is of weight 2, the Jacobian J 1 ( 16) of X 1 (16) contains the mod l representation ρ f,l attached to l for any prime l of K f .For this example, let us take l = (5, i − 2), one of the two primes of K f above 5.
As explained in Sect.3.3, we must begin by choosing a prime p to work with.After trying all primes up to 100, which requires computing a p (f ) for p 100 by classical methods, we decide to take p = 23, because ρ f,l (Frob 23 ) has order only 4.
Next, we construct a 23-adic Makdisi model of X 1 ( 16) with residue degree a = 4 and accuracy O(23 e ), where we have chosen e = 7 so as to identify rationals of height at most 4 × 10 4 .This construction involves spotting the elliptic curve which has all its 16-torsion defined over the degree-4 unramified extension of Q 23 .
We then generate an F l -basis of the subspace of J 1 ( 16)(F 23 ) [5] that affords ρ f,l .On our way, we confirm that the rational canonical form of ρ f,l (Frob 23 ) is 0 −2 ∈ GL 2 (F l ); this was the only possibility, since during the first step, we had determined from the value of a 23 (f ) mod l that the characteristic polynomial of ρ f,l (Frob 23 ) is x 2 + 2x + 2, which is separable mod 5.
We must now lift a basis of the representation space to J 1 (16)(Z 23 4 /23 7 ) [5].Actually, since we know that the action of Frob 23 on the representation space is cyclic, we can afford to only lift one 5-torsion point.We then generate all the points of the representation space over Z 23 4 /23 7 by mixing the group law of the Jacobian and the action of Frob 23 , and we evaluate the resulting points by 20 versions of the evaluation map α, compute the corresponding 20 versions of the polynomial F (x), and keep the nicest one.We thus find that our Galois representation is described by the polynomial Remark 6.1.1Since the computation also returns an indexation of the 23-adic roots of F (x) by the nonzero vectors of F 2 l , we can easily compute a polynomial describing the projective version if we wish to do so, by gathering symmetrically (in this instance, summing) the roots of F (x) along the vector lines of F 2 l .We find the polynomial which has one rational root and one irreducible factor of the degree 5.The representation ρ f,l is thus reducible, a fact that can easily be checked independently.

mod 13
As a second example, we compute the representation attached to = 1.12.a.a = q − 24q 2 + 252q 3 + O(q 4 ) mod l = 13.This representation ρ ,13 is found in the 13-torsion of the Jacobian J 1 (13) of the modular curve X 1 (13), whose genus is again 2. Nevertheless, its computation requires significantly more effort than the previous case, because we have to have to process 13 2 − 1 = 168 torsion points and thus output a polynomial F (x) of degree 168, as opposed to 5 2 − 1 = 24 in the previous example.Besides, whereas in the previous example the image of the presentation was contained in a Borel subgroup of GL 2 (F 5 ), the image of ρ ,13 is the whole of GL 2 (F 13 ), which is considerably larger; as a result, we expect it to be more difficult to find a prime p such that ρ ,13 has small order, so this time we try primes up to 200.However, this turns out not to be necessary: indeed, for p = 73, the order of the image of the Frobenius is again a = 4 only.
As deg F (x) = 168 and as we expect its coefficients to be less nice than in the previous example, we choose to work at accuracy O(73 44 ), so as to identify rationals of height up to 10 40 .We construct a 73-adic Makdisi model of X 1 (13) with residue degree a = 4 at this accuracy by spotting the elliptic curve which has all its 13-torsion defined over the degree-4 unramified extension of Q 73 .
We then generate an F 13 -basis of the subspace of J 1 (13)(F 73 ) [13] that affords ρ ,13 , and we lift to accuracy O(73 44 ) a 13-torsion point which generates the representation space under Frob 73 .Finally, we generate all the points of the representation space over Z 73 4 /73 44 , evaluate them by several versions of the map α ∈ Q J 1 (13) , and keep the nicest resulting polynomial F (x).
In the end, we find that our Galois representation is described by a polynomial of the form whose coefficients have 10103 2 as a common denominator, and numerators of up to nearly 40 decimal digits.

mod 19
We now compute the representation attached to f = again, but this time mod l = 19.This involves processing 19 2 − 1 = 360 torsion points, but more importantly, this computation takes place in the Jacobian of a curve of genus g = 7, namely X 1 (19), and is therefore significantly more difficult than the previous ones.
After having tried all primes p 1000, we select p = 107, since it allows us to work in residue degree a = 6.We thus construct a 107-adic Makdisi model of X 1 (19) with degree 6, and accuracy O (19 247 ) so as to identify rationals of height up to 10 250 , and obtain a polynomial F (x) of degree 360 describing ρ ,19 .

Timings and comparison with the complex method
The table below summarises the parameters and the execution real times (as opposed to CPU times) on the author's laptop of the three computations described above.We specify these times by decomposing the computation into five phases: • The choice of a prime p, by trying all primes up to some bound and keeping the value for which the order of the image of the Frobenius at p is minimal, • The construction of a p-adic Makdisi model of the modular Jacobian J , • The determination of an F -basis of the representation space in J (F q )[ ], • The p-adic lift of this basis from J (F q ) to J (Z q /p e ), • The generation of all the points of the representation space in J (Z q /p e ) and the evaluation of (several versions of) a rational map α ∈ Q(J ) at these points.As a comparison, a few years ago, the computation [3] of ρ ,13 and of ρ ,19 by the complex-analytic method took about 5 minutes and 40 minutes respectively, even though they took place on 64 cores on the supercomputing cluster [24] instead of the 4 cores of the author's laptop.
Our implementation of the p-adic method thus strongly outperforms our implementation of the complex-analytic method [3], and makes it reasonable to compute modular Galois representations on a personal laptop (even though having many cores would have sped up significantly the computation of ρ ,19 as it involves processing 360 torsion points).That we implemented our p-adic approach in C language whereas [3] was written in Python probably plays a part in this, but there are other, more fundamental reasons for this difference of performance.
Indeed, in order to generate -torsion points, the complex-analytic method begins by computing a high-accuracy approximation over C of a period lattice of the modular curve, which takes a significant amount of time since it requires, in particular, computing the q-expansion to high accuracy of a basis of the space of cusp forms of weight 2. On the contrary, the p-adic lifting method starts in "low accuracy", meaning mod p, where torsion points can be obtained easily thanks to fast exponentiation; therefore it does not suffer from this overhead.
Besides, as explained in [7, 6.4], since the evaluation map α from the Jacobian to A 1 is by design Galois-equivariant, the p-adic method can save a lot of effort by computing and evaluating not all the points of the representation space, but only one per orbit under the Frobenius Frob p .In contrast, the complex method can only use complex conjugation, which has order 2 and thus can only halve the amount of work.
Another pleasant feature of the p-adic approach is that it can naturally deal with nonsquarefree levels, as demonstrated by the first example above which took place in level N = 16.On the contrary, non-squarefree levels are problematic for the complex method, since the computation of the periods of the modular curve requires the determination of Atkin-Lehner pseudo-eigenvalues [4, 3.2.3],which cannot in general be easily read off the coefficients of a newform when the level is not squarefree [4, 2.1.2].Remark 6.2.1 Experimenting shows that if we had used evaluation maps α from the Jacobian to A 1 based on rational forms instead of q-expansions (see Sect. 5), we would have had to increase the p-adic accuracy from O(23 7 ) to about O (23 600 ) in the first example, and from O(73 44 ) to O(73 8000 ) in the second example, which would have made the computations tremendously slower.

"Larger" examples
To conclude, we demonstrate the performance of our method on "larger" computations by providing in the table below the execution real times of three examples on the supercomputing cluster [24].The first of these representations is attached to a newform which admits a companion form, and the second one is attached to a supersingular newform, which is the reason why we had already computed these representations by the complexanalytic method [3] in [4].We had computed the third representation in [3]  These timings show that our p-adic approach remains suitable for "large" computations.In particular, we observe that the construction of a high accuracy p-adic Makdisi model of the modular Jacobian is very far from being the bottleneck of the computation of the representation, and that the p-adic lifting method which we introduced in [7] remains efficient in high genera.In comparison, the computation of these representations by the complex-analytic method [3] had respectively taken 12 hours, 12 hours, and 5 days.that commutes for all i.Then, by the same line of ideas as for Frob p , we may apply ϕ * to a matrix W D x representing a point x ∈ J : all we have to do is permute its rows by σ −1 φ .Again, this takes negligible time compared to group operations in J .
Indeed, on the one hand, as x ∈ J is represented by L(−D x ), its image ϕ * (x) is represented by L(−ϕ * (D x )) since L is invariant by ϕ; and on the other hand, each column of W D x consists of the vector of values (under the t i ) of a global section w of L ⊗2 (−D x ), so permuting its entries by σ −1 ϕ yields the vector of values of a section w such that w (P σ (i) ) = w(P i ), so that indeed w = N ϕ (w).
This idea allows us to implement the action of the diamond operator d on Makdisi models of modular curves X H (N ).

A.2 Fast exponentiation using cyclotomic polynomials and the Frobenius
Let again C be a curve over Q, and suppose that the -torsion of its Jacobian J contains a Galois representation that we wish to compute.
In the notation of Strategy 3.1.10,we typically have M ≈ #J (F q ) ≈ q g = p ag .Therefore, for large genus g, on the one hand M is quite large, especially as (3.1.8)usually imposes that q is in the thousands or even millions; and on the other hand, performing one addition in J using Makdisi's algorithms relies on linear algebra of size O(g) over F q , which is rather costly.Thus, even though we use fast exponentiation, multiplication by M in J (F q ), which is required to generate -torsion points, can take a significant amount of time.
However, we have seen in the previous Subsection that applying the Frobenius Frob p takes negligible time, so it is natural to try to use the action of Frobenius to speed up this multiplication-by-M step.For this, we begin by establishing the following result: Lemma A.2.1 Let G be a finite Abelian group, and let φ : G → G be an endomorphism.View G as a Z[x]-module with x acting as φ.Suppose we know a monic polynomial P(x) ∈ Z[x] such that P(φ) = 0 and which factors in Z[x] as P = AB with Res(A, B) coprime to #G (in particular, A and B must be coprime in Q are inverses of each other. Suppose now that a; if one has chosen a minimal for the points of T to be defined over F q = F p a , which is what we will do in practice, then this is equivalent to saying that ρ T (Frob p ), which has order a, is semisimple.Take G to be the part of J (F q ) coprime to a (which has thus the same -torsion as J (F q )), φ = Frob p : G → G, and P(x) = x a − 1, which factors over Z into the cyclotomic polynomials Since #G is prime to disc(x a − 1) = ±a a by construction, an iterated use of Lemma A.2.1 yields the decomposition so each of these factors is typically much smaller than J (F q ).
We can thus obtain -torsion points as follows: 1. Pick d | a such that | N d (there will be at least one).

Write
3. Take a random x ∈ J(F q ).Apply x a −1 Φ d (x) (Frob p ) to it, multiply the result by M d , and then repeatedly by until we get 0. Return the last nonzero point.

Strategy A.2.2 Generating torsion points thanks to cyclotomic polynomials
This method is valid since as a by assumption, M d is divisible by [J (F q ) : G], so multiplication by M d includes the effect of projecting from J (F q ) to G. Its advantage is that the number of required operations in J (F q ) is approximately lg N d ≈ ϕ(d)g lg p, compared with lg N ≈ ag lg p with Strategy (3.1.10).Indeed, typically the cofactor x a −1 d (x) has few nonzero coefficients, and these coefficients are usually ±1, so applying x a −1 d (x) (Frob p ) requires few operations in J (F q ) and thus takes negligible time since applying Frob p is effortless.
Example A.2.3 Let C be the modular curve X 1 (13), which has genus 2, and let J be its Jacobian.Suppose we want to generate -torsion points of J where = 29, for example.Using formula (3.2.1) below and [7, Proposition 5.1], we find that for p = 191, J [ ] is defined over F p a for a = 12.We have where lg M ≈ 162, so if we use the method presented in Strategy (3.1.10),then we need to perform about 200 additions in J (F p a ) to obtain an -torsion point.
In comparison, if we take d = 12, we find where lg M 12 ≈ 60, so we can produce an -torsion point with fewer than 100 additions in J (F p a ) by using Strategy (A.2.2), even taking into account the operations required to multiply by Similarly, for d = 3 we have so this is the only d the rest of J [ ] comes from, and we have lg M 3 ≈ 20 only.However, this produces -torsion points defined over F p 3 , so if we want to get all of J [ ], then we need to generate points using d = 12 as well.
Remark A.2.4 In many cases, we do not only want points of J [ ], but actually points in the piece T χ of J [ ] where Frob p acts with some prescribed characteristic polynomial χ(x).Therefore, in Strategy (A.2.2), we should only consider the d | a such that d (x) has a nontrivial common factor mod with χ(x).

A.3 Optimising the Makdisi model
Makdisi's algorithms to compute in Jacobians rely on linear algebra involving matrices whose dimensions are determined by the parameters d 0 and n Z introduced in Definition 3.1.1.The smaller these parameters, the faster the computations; however these parameters must respectively satisfy the bounds (3.1.2) and (3.1.3)for Makdisi's algorithms to be valid.The p-adic Makdisi model of X H (N ) that we constructed in Sect.3.2 often significantly exceeds these bounds.We are now going to show how to optimize it, by tweaking it so that (3.1.2) and (3.1.3)are satisfied as sharply as possible.In practice, this optimisation results in a major speedup.
Let us begin with (3.1.2).Recall that we chose L to be the line bundle whose sections This ensured that the bound (3.1.2) on d 0 is satisfied since we assume that the number ν ∞ of cusps of X H (N ) is at least 3.We can thus make (3.1.2) an equality as in [3, 3.3], that is to say, by fixing three cusps of X H (N ) and by replacing L with the sub-bundle whose section space consists of the modular forms of weight 2 that vanish at all cusps except maybe these three.
This space may be determined by linear algebra, provided that we can determine the value of modular forms ranging over a basis of M 2 H (N ) at every cusp.Note that these values are well-defined by (2.3.13)applied to the case n = 0.
Recall that in Sect.4.4, we found a basis of M 2 H (N ) consisting of forms f v i ,w i 2,H as defined by (4.4.1) for various v i , w i ∈ (Z/N Z) 2 .Since f v 1 γ = f vγ 1 for all γ ∈ SL 2 (Z/N Z), it is actually enough to determine the value of f v 1 at the cusp ∞ in terms of v, which we can do thanks to Theorem 4.2.2.
We can thus optimise L so that the bound (3.1.2) on d 0 is sharp.Since d 0 = deg L drops, we can then also reduce the number of points in the fiber of X H (N ) −→ X(1) at which we evaluate our forms.In fact, by (3.1.3),we only need to retain 10g + 6 of these points, since we now have d 0 = 2g + 1 exactly.However, the Definition 3.1.1 of a Makdisi model also requires the set of these points to be globally invariant under Frob p .Since we have determined how Frob p acts on the N -torsion of the elliptic curve corresponding to this fibre, we can explicitly decompose this fibre into Frob p -orbits; we then discard some of these orbits so that the bound (3.1.3)is satisfied and as sharp as possible, and we only evaluate our forms at the points in the remaining orbits to construct our p-adic Makdisi model of X H (N ).We can usually make (3.1.3)almost sharp since we have chosen the prime p so that the order a of Frob p is small.

A.4 Finding a suitable elliptic curve and efficiently computing its N-torsion
Let N ∈ N be an integer, p 6N a prime, and a ∈ N a degree such that there exists an elliptic curve E defined over Q, having good reduction at p, j(E) ≡ 0, 1728 mod p, and E[N ] defined over the unramified extension Q q of Q p of degree a; in particular, the bounds (3.3.1) and (3.3.2) must be satisfied.We are going to present efficient methods to find such a curve efficiently, and to compute data needed to derive a p-adic Makdisi model of modular curves, such as high-accuracy p-adic approximations the coordinates of E[N ].The typical range that we have in mind is N 1000, p 10 4 , and a 100.
Since p N , the requirement that ) ⊂ F q , where Ē denotes the reduction of E mod p; therefore, since p 6, we will actually look for integers 0 < A, B < p such that the short Weierstrass equation (W) : y 2 = x 3 + Ax + B viewed mod p defines an elliptic curve ĒA,B over F p (so 4A 3 + 27B 2 ≡ 0 mod p).Since A, B ≡ 0 mod p, j( ĒA,B ) / ∈ {0, 1728}.Our assumption is that there exists at least one such pair (A, B) such that ĒA,B [N ] is defined over F q , so that (W) then defines an elliptic curve over Q with the desired properties.
Our strategy simply consists in trying random pairs (A, B) until the condition is satisfied.Given such a random pair, we expect that (A.4.1) will most likely not be satisfied, so instead of directly testing (A. for all k.If we let Frob p : x → x p be the standard pro-generator of Gal(F p /F p ), and if we define Frob q = Frob a p : x → x q , then this can be rephrased by saying that Frob q must act trivially on ĒA,B [l Since by assumption p is not too large, given a pair (A, B), we can quickly determine the quantity The characteristic polynomial of Frob p acting on ĒA,B is then Given a p , it is thus straightforward to compute its discriminant = a 2 p − 4p ∈ Z, as well as the Newton sum where α and β are the roots of χ(x) in Q and a = [F q : F p ] as above; this can even be done symbolically, without actually computing α and β.Naturally, if the cardinality # ĒA,B (F q ) = q + 1 − ν a is not a multiple of N 2 , then the pair (A, B) must be rejected. Define then the action of Frob q on ĒA,B [M 1 ] is semisimple and has characteristic polynomial else, Frob q ĒA,B [l k ] is unipotent iff a a p ≡ 2 a mod l k , therefore the pair (A, B) can be rejected if this condition is not satisfied for at least one of such l k .
These three simple tests eliminate most of the pairs (A, B).We now assume that (A, B) has passed these three tests, so the action of Frob q is trivial on ĒA,B [M 1 ] and unipotent on ĒA,B [l k ] for each l k | ; it remains to determine whether Frob q really acts trivially on ĒA,B [l v k k ] for each k.This is automatically the case for the l k such that v k = 1, as well as for the l k | such that v k = 1 and l k | a since a unipotent mod l k matrix of size 2 × 2 which is also an a-th power is then necessarily trivial; we therefore do not consider these primes anymore.
For each of the remaining primes, we then compute the division polynomial of ĒA,B , and determine the degrees of its factors over F p , which is faster than factoring it completely [26, 3.4.3].If these degrees do not all divide a, then this polynomial does not split over F q , so the action of Frob q on ĒA,B [l v k k ]/ ± 1 is nontrivial and the pair (A, B) can be rejected.Else, for each k such that l v k k = 2, we determine the roots of ψ l v k k (x) in F q , and for each such root z, we check whether z 3 + Az + B is a square in F q by raising it to the q−1 2 using fast exponentiation (if l v k k = 2, then this will automatically be satisfied since ĒA,B [2]/ ± 1 = ĒA,B [2]).If this is the case, we have found a suitable pair (A, B).
Suppose now that we have found a suitable pair (A, B).In order to compute a p-adic Makdisi model for X H (N ) to accuracy O(p e ), where e ∈ N is a fixed parameter, we need to determine the coordinates in Z q /p e of the N -torsion points of the elliptic curve E A,B over Z p defined by (W).It is sufficient to determine the coordinates of two points P, Q forming a basis of E A,B [N ], since the coordinates of the other torsion points can then be obtained by applying the group law of E A,B (Z q /p e ) as (W) has good reduction at p. We must also compute the matrix expressing how Frob p acts on E A,B [N ] with respect to this basis, so as to determine how Frob p permutes the points of the fibre of X H (N ) −→ X (1) corresponding to E. Finally, we need the value e N (P, Q) of the Weil pairing of this basis so as to use the q-expansion formula established in Theorem 4.2.2.
Again, we want to try to avoid the expensive computation of the N -division polynomial of E A,B , so we proceed prime-by-prime.Factor as above N = k l v k k , and let i k ∈ Z/N Z be the idempotents corresponding to the Chinese remainder decomposition For each k, we begin by computing the polynomials , where ȳP k ∈ F q is either square root of x3 P k + Ax P k + B, and similarly for ȳQ k .Then Pk , Qk ∈ ĒA,B (F q ) are two points of exact order l k ] over F q .We assume that this is the case from now on.
We can then determine the matrix of Frob p acting on E A,B [l v k k ] with respect to (the unique p-adic lift of) this basis as where log zk : k Z denotes the discrete logarithm in base zk .Next, we lift this basis ( Pk , Qk ) from F q to Z q /p e by first Hensel-lifting the x-coordinates as roots of , and then the y-coordinates as square roots of x 3 + Ax + B; we thus obtain a basis In principle, the value z k ∈ Z q /p e of its Weil pairing could also obtained by lifting zk as a root of x l v k k − 1, but we defer this for now since we will see that we can do better.
Remark A.4.2 Some of the division polynomials ψ l v k k (x) may have been computed in F p [x] during the earlier phase when we searched for an appropriate pair (A, B); however, they need to be re-computed, since we need their value in Q[x] (as opposed to mod p) here.

It is then clear that
, with respect to which the matrix ∈ GL 2 (Z/N Z) describing the action of Frob p can be obtained from the k by Chinese remainders thanks to the idempotents i k .Besides, its Weil pairing z = e N (P, Q) ∈ μ N (Z q /p e ) (A.4.3) may be determined off the zk .Indeed, it is enough to determine z = z mod p ∈ F q , since we can then Hensel-lift this value as a root of x N −1.Furthermore, given an elliptic curve E and two integers M 1 , M 2 ∈ N, the definition of the Weil pairing in terms of meromorphic functions on E with prescribed divisors shows that we have the identity The advantage of this approach is that the Weil pairing computations, which can be expensive, are only performed on the ĒA,B [l

A.5 Using the Hecke operator T p A.5.1 Necessity of the use of T p
Let f = q+ n 2 a n (f )q n ∈ N k (N, ε f ) be a newform, and suppose that we wish to compute the Galois representation ρ f,l attached to f mod a prime l of degree of ρ f,l (Frob p ) is coprime mod with its cofactor in the local L factor In some rare cases, it may happen that condition (3.1.9)is not satisfied by any prime p, so our p-adic method to compute ρ f,l , as presented so far, does not apply.
Example A.5.1 Let f be the newform of [23] label 5.6.a.a and l = 13.Then χ p (x) 2 | L p (x) mod for all p = 13.This phenomenon is related to f being supersingular at l and having trivial nebentypus.
However, in the immense majority of cases, multiplicity one statements such as [9,Theorem 3.5] show that this can be remedied by isolating the subspace T f,l of dimension 2 of J H (N )[ ] affording ρ f,l as the subspace where the Hecke algebra acts with the eigenvalue system of f mod l instead of in terms of the action of the Frobenius at a good prime p.In other words, the representation space T f,l can be carved out as yield an F -subspace of dimension 4 which is a non-split extension of one copy of ρ f,l by another.This explains why an attempt based on the characteristic polynomials χ p (x) alone cannot succeed in this case.
In order to remedy this situation, (A.5.2) thus suggests we implement the action of the Hecke algebra on Makdisi models of modular curves.Although this is theoretically possible since [19] shows that pull-backs and push-forwards are computable in Makdisi models, this approach seems complicated, and we have chosen not to follow it.Instead, we hope that there exists a prime p such that the representation space can be carved out simply as (A. 5.3) which is tantamount to having a p (f ) mod l = a p (f ) mod l for a unique pair (f , l ) with f an eigenform of weight 2 and level H (N ) and l a prime of K f above .We have not encountered any case where no such p exists, so this seems to be a reasonable assumption.

A.5.2 Implementing T p on p-adic Makdisi models
The Eichler-Shimura relation [13, 8.7.2] states that Therefore, if we construct a p-adic Makdisi model of J H (N ) with the same prime p as in (A.5.3) and which contains the extra data required to apply p * , then we get an implementation of T p on J H (N )(F p ), and we may thus alter Strategy 3.1.10to carve out the representation space with T p instead of Frob p .In order to construct such a p-adic Makdisi model, we must first find a good prime p such that (A.5.3)holds.We follow a search procedure similar to that described in Sect.3.3, where we substitute (A.5.3) for the condition that χ p be coprime with its cofactor in L p .In order to test (A.5.3) for a given prime p, we could compute the matrix of T p (with respect to any Z-basis) acting on cuspidal modular symbols of level H (N ), reduce it mod , and see if its a p (f ) mod l-eigenspace has dimension 2 and not more.Alternatively, in view of the pairing between cusp forms and cuspidal modular symbols, we can also compute the matrix of T p acting on the space S 2 H (N ) of cusp forms with respect to a basis of S 2 H (N ) such that this matrix has integers entries, which can be done with [11], reduce it mod , and see if its a p (f ) mod l-eigenspace has dimension 2 and not more.Indeed, although it is conceivable that the basis of S 2 H (N ) chosen by [11] does not remain a basis after reduction mod , or that the pairing between modular forms and modular symbols degenerates mod , these phenomena can only increase the dimension of the a p (f ) mod l-eigenspace, so we do get a sufficient criterion for (A.5.3) to be satisfied.The advantage of working with cusp forms rather than modular symbols is that we will need anyway to determine the matrix of T p acting on cusp forms when we compute L p (x) by (3.2.1).
Next, our p-adic Makdisi model must satisfy the extra requirements outlined at the end of Sect.A.1 so we can apply the diamond operator p * .The permutation induced by p on the fibre of X H (N ) −→ X( 1) is easily determined thanks to (2.1.5).However, we must modify the optimisation process described in Sect.A.3, for two reasons.First, the line bundle L must be invariant by p , so instead of taking the bundle whose sections are the forms of weight 2 that vanish at all but three cusps, we must take the bundle whose sections are the forms of weight 2 that vanish at all cusps outside a set S containing at least three cusps and which is stable not only by Gal(Q/Q), but also by p .Therefore, we may need to take a set S with slightly more than three elements, which makes the d 0 = deg L larger and thus takes us away from attaining sharpness in the bound (3.1.2).Second, while we can still drop some of the points of the fibre of X H (N ) −→ X(1), we must ensure not only that the bound (3.1.3)is satisfied in spite of the accretion of d 0 , but also that the remaining points are globally invariant not only under Frob p , but also under p .As a result, we may not be able to optimise our p-adic Makdisi model as sharply as before, in which case all our computations in J H (N ) will be a little slower.For this reason, it is preferable to carve out T f,l with Frob p as before when possible, and to reserve this new approach to tenacious cases such as Example A.5.1.
Remark A. 5.3 In principle, it would be possible to construct two p-adic Makdisi models of X H (N ), namely a "large one" with the extra data for p * and a "small", better-optimised one without.We would then use the large model to generate points of T f,l , convert them to points of the small model, and proceed with the small model for the rest of the computation.We have not yet implemented this feature.

A.5.3 Isolating T f,l by T p
We now describe in detail our new approach to compute ρ f,l , given a good prime p satisfying (A.5.3).Although this approach remains valid even if χ p (x) L p (x), we are chiefly interested in the case where χ p (x) ∦ L p (x).
As before we can determine from χ p (x) an integer a ∈ N such that the points of T f,l and the -th roots of unity are defined over F q = F p a .We then construct a p-adic Makdisi model of J H (N ) over Q q which contains the extra data needed to apply p * .We then reduce this model mod p to obtain a Makdisi model of J H (N )(F q ), while retaining the high p-adic accuracy model for later use.Writing J = J H (N )(F q ) for brevity from now on, we can then apply T p on J as where m is any integer satisfying m ≡ p mod .In particular, working with large values of p does not slow down the application of T p , so there is no harm in choosing a large value of p so as to make a smaller, thus making the calculations in J faster.
Let υ p (x) be the largest mod factor of L p (x) whose irreducible factors all divide χ p (x); by construction, χ p (x) | υ p (x) L p (x), so we know how to generate random points of the subspace U = ker υ p (Frob p ) ⊆ J [ ] by using the action of Frob p .Observe that U ⊇ T f,l is an F -space of dimension d = deg υ p (x) since υ p (x) L p (x).
As F q contains the -th roots of unity by assumption, Algorithm 6 of [7] allows us to evaluate the Frey-Rück pairing where the rightmost arrow consists in x −→ x (q−1)/ followed by the discrete logarithm with respect to some fixed primitive -th root of unity.This paring is perfect, so we will use it to detect F -linear dependency in J [ ], similarly to Algorithm 13 of [7].Consider now the following algorithm.Then set x ← T p x, and go to step 3.
5. Else, let λ 1 , • • • , λ r+1 ∈ F be the coefficients of a non-trivial linear dependency relation between the columns of P .
6.If r i=1 λ i x i +λ r+1 x does not evaluate to 0 ∈ J, then first set d + ← d + +1, r ← r + 1, x r ← x, and P ← P .Next, generate random points y ∈ J until the row [x 1 , y], • • • , [x r , y] is linearly independent from the rows of P , find an i d such that the i-th row of P is in the F -span of the other rows of P , and replace this i-th row with [x 1 , y], • • • , [x r , y] and y i with y.Finally, set x ← T p x, and go to step 3.
8. If r i=1 λ i x i + λ r+1 x = 0 ∈ J but d + = 0, set λ i ← −λ i /λ r+1 for all i r so that x = r i=1 λ i x i .Replace M with ⎛ where the block with the sub-diagonal of ones has d + rows, and then determine K = ker(M − a p (f ) mod l).If dim K = 2, then let (κ 1,j , • • • , κ r,j ) (j ∈ {1, 2}) be a pair of vectors of F r forming a basis of K, let b j = r i=1 κ i,j x i , and return the pair (b 1 , b 2 ) of points of J[ ].Else, go to step 2.

Strategy A.5.7 Carving out T f,l by T p
Every time we enter step 2, the r < d points x 1 , • • • , x r ∈ U are linearly independent over F and span a subspace X ⊆ U which is stable under T p , M is the matrix of T p on X with respect to x 1 , • • • , x r , and the d × r matrix P of pairings [x j , y i ] has full rank r.Then at step 3, the d × (r + 1) matrix P has rank either r + 1 or r.If P has rank r + 1, then x is linearly independent from the x j , so at step 4 we append x to the x j , increase r, and start over with T p x instead of x.If P has rank r, then at step 5 the λ i are unique up to scaling, and either x is genuinely linearly dependent on the x j , or the linear forms [•, y i ] do not separate the points of X; we remove this ambiguity by checking whether r i=1 λ i x i +λ r+1 x is genuinely zero.If not, then at step 6, we append x to the x j thus increasing r, then we modify one of the y i so that the linear forms [•, y i ] separate the points in the span X of the x j , and finally we start over with T p x instead of x.If r i=1 λ i x i + λ r+1 x = 0, then x, is linearly dependent on x 1 , • • • , x r since the latter are linearly independent; in particular, λ r+1 = 0.At this stage, x 1 , • • • , x r span a T p -stable F -subspace X of U , whose dimension r has increased by d + since the last time we entered step 2. So if d + = 0, then the point x generated at step 2 was already in X at that time, so we simply try again with another random x ∈ U .Else, at step 8 we update M, and see whether X contains the 2-dimensional subspace T f,l = ker(T p − a p (f ) mod l).If it does, then we return a basis of T f,l ; else we go back to step 2 to enlarge the F [T p ]-module X by throwing in a new random point x ∈ U .
We thus obtain a pair of points (b 1 , b 2 ) of J [ ] forming an F -basis of the representation space T f,l , so we may proceed with the calculation of ρ f,l by our usual Strategy 3.1.10from step 4 on.

A.5.4 An explicit example
We computed in Sect.6.3 the mod = 13 representation attached to the newform f = 5.6.a.a.We observed in Example A.5.1 that the piece of the modular Jacobian J which affords it cannot be isolated by the action of Frob p for any prime p, so we were forced to use the variant relying on the action of T p .
A search over the primes p 1000 with the methods presented in Sect.A.5.2 reveals that the representation space may be isolated as T f,13 = ker T p − a p (f )| J [13] (A.5.8) for many primes p.We choose p = 439, which satisfies (A.5.8) and has the extra advantage that the order of the image of Frob p is a = 6 only.
Having constructed a p-adic Makdisi model of J as usual, our implementation then generates a few random points of J [13](F 439 6 ) by the method outlined in Sect.A.2.The first of these points, which was generated in the subspace killed by 3 (Frob 439 ), spans an F 13 [T 439 ]-module S 1 of dimension dim  (2.3.10), and yet the matrix M = 1 2 2 5 ∈ SL 2 (Z/N Z) takes ∞ to s and yields rational q-expansions by (2.3.17).
One also easily verifies that this converse also fails for all the values of N such that φ(N ) < 3, apart from N = 1 for which it is trivially true.

1 . 6 )
. In particular, two cusps (c, d), (c , d ) of X 1 (N ) are in the same Galois orbit iff c = ±c mod N .One easily verifies that the correspondence between this representation of the cusps by pairs (c, d) and the more traditional one by classes of elements of Q ∪ {∞} is as follows: given an element a/c ∈ Q ∪ {∞} in lowest terms, find b, d ∈ Z such that a b c d ∈ SL 2 (Z); then the cusp represented by a/c in the traditional notation is represented by (c, d) in our representation, and vice-versa.

Remark 3 . 1 . 4
Makdisi's algorithms deal with global sections of powers L ⊗n with n up to 5. The bound (3.1.3)ensures that such sections are faithfully represented by their values at the points P i .

Theorem 4 . 2 . 2
Let z = e N β(1, 0), β(0, 1) be the primitive N -th root of unity which indicates the geometric component of X(N ) we are working in.There exists a constant C depending only on N such that for all 0 c < N and d

4 . 1 )
by computing the N -division polynomial ψ N (x) of ĒA,B which would be time-consuming, we begin by submitting ĒA,B to a battery of quick tests based on point-counting and aiming at weeding out most of the pairs (A, B) for which (A.4.1)does not hold.Once we find a pair (A, B) which passes these tests, we then submit it to extra tests to try to prove that (A.4.1)does hold, still while trying to avoid expensive computations such as the determination of ψ N (x).Since N is reasonably small, we can factor it as N = k l v k k where the l k ∈ N are distinct primes, and then (A.4.1) is equivalent to

= 1 .
If this is not the case, then we start over with another choice of xP k , xQ k ; else we have obtained a basis of ĒA,B [l v k

T p =
Frob p +p p * Frob a−1 p .(A.5.5)Remark A.5.6 We will actually need to apply T p to points of J [ ] only.Therefore, we can replace (A.5.5) with T p = Frob p +m p * Frob a−1 p

1 .
Generate random points y 1 , • • • , y d of J[ ].Initialise a matrix M of size 0 × 0 over F , and a matrix P of size d × 0 over F .Finally, set r ← 0. 2. Set d + ← 0, and generate a random point x of U . 3. Set P ← the matrix obtained by sticking the column containing [x, y 1 ], • • • , [x, y d ] to the right of P .4. If P has rank r + 1, then set d + ← d + + 1, r ← r + 1, x r ← x, and P ← P .

F0 1 8 .
13 S 1 = 2, on which the matrix of T 439 is 0 Since a 439 (f ) = 8 mod 13, this module S1 does not yet contain the 2-dimensional space T f,13 = ker(T 439 − 8 mod 13), so we enlarge S 1 by including the second random 13-torsion point, which was generated in the subspace killed by 2 (Frob 439 ).The dimension of the resulting F 13 [T 439 ]-module Finally, if c ≡ 0, then 1 = gcd(d, N ), and (A.5.1)implies that for all x ∈ (Z/N Z) × , d(x+1) or d(x − 1) vanishes mod N .Since d is invertible mod N , this implies x ≡ ±1, so again this case cannot occur.Corollary A.6.4 Let s be a cusp of X 1 (N ).If s has width N , then there exists M ∈ SL 2 (Z/N Z) such that M • ∞ = s and that M yields rational q-expansions.The converse holds, provided that φ(N ) 3 and that N ≡ 2 mod 4.ProofLet (c, d) represent the cusp s, and suppose s has width N .Then gcd(c, N ) = 1 by Proposition A.6.1.Let (c, d) represent the cusp s.Since d lives in (Z/(c, N )Z) × , we may replace it with d t = d + tc for any t ∈ Z.As gcd(c, N ) = 1, we can choose t so that d t ≡ 0 mod N , and then (2.3.17)obviously holds.Conversely, suppose there is such an M = a b c d ∈ SL 2 (Z/N Z), so that s is represented by (c, d), and that N | 2d by Proposition A.6.2.We are going to deduce that gcd(c, N ) = 1, which will conclude by Proposition A.6.1.Observe that 1 = gcd(c, d, N ).We distinguish two cases: if N is odd, then N | d, so gcd(c, N ) = 1; and if 4 | N , then N 2 | d, so gcd(c, N /2) divides gcd(c, d, N ) = 1 and is therefore 1, and since N /2 is even, this implies gcd(c, N ) = 1.Remark A.5.2 The converse part of Corollary A.6.4 may fail if N ≡ 2 mod 4, as illustrated by the example N = 10, s represented by (c = 2, d = 5).Indeed, this cusp has width w = 5 < N by acw , so we want c 2 w ≡ 0 mod N and 1±acw ∈ H (note that (1+acw)(1−acw) = 1−a 2 c 2 w 2 ≡ 1 if c 2 w ≡ 0, so the ± sign is irrelevant, in that this identity with either sign implies the one with the other sign).Let g 2 = gcd(N, c 2 ) and N 2 = N /g 2 ; then c 2 w ≡ 0 iff N 2 | w, whence finally w = N 2 min{t ∈ N | 1 + acN 2 t ∈ H}.
3.2, we had left open the problem of evaluating (with respect to some appropriate local trivialisation) a basis of the space M 2 H (N ) of global sections of the line bundle L of our p-adic Makdisi model of J H (N ) at points of X H (N ) represented as pairs (E, H • Q), where E is a fixed elliptic curve and Q varies in E[N ].We now show how this can be done thanks to Makdisi's forms f v 1 .For simplicity, we will arrange that E is specified by a Weierstrass equation (W) : The computation from the values of α of the polynomial F (x) which describes the representation takes negligible time.
1 mod N by assumption, so this action is trivial iff ν a ≡ 2 mod M 1 .We can therefore reject the pair (A, B) if this condition is not satisfied.If now l k is one of the prime factors of N dividing , then Frob p ĒA,B [l k ] is not semisimple, and therefore has a single eigenvalue c ∈ F l k which satisfies 2c ≡ a p mod l k .If l k = 2, then Frob p ĒA,B [l k ] is necessarily unipotent, and therefore so is Frob q ;