
// Make newform 3.10.a.b in Magma, downloaded from the LMFDB on 25 April 2024.
// To make the character of type GrpDrchElt, type "MakeCharacter_3_a();"
// To make the coeffs of the qexp of the newform in the Hecke field type "qexpCoeffs();"
// To make the newform (type ModFrm), type "MakeNewformModFrm_3_10_a_b();".
// This may take a long time!  To see verbose output, uncomment the SetVerbose lines below.
// The precision argument determines an initial guess on how many Fourier coefficients to use.
// This guess is increased enough to uniquely determine the newform.
// To make the Hecke irreducible modular symbols subspace (type ModSym)
// containing the newform, type "MakeNewformModSym_3_10_a_b();".
// This may take a long time!  To see verbose output, uncomment the SetVerbose line below.
// The default sign is -1.  You can change this with the optional parameter "sign".
function ConvertToHeckeField(input: pass_field := false, Kf := [])
    if not pass_field then
        Kf := Rationals();
    end if;
    return [Kf!elt[1] : elt in input];
end function;


// To make the character of type GrpDrchElt, type "MakeCharacter_3_a();"
function MakeCharacter_3_a()
    N := 3;
    order := 1;
    char_gens := [2];
    v := [1];
    // chi(gens[i]) = zeta^v[i]
    assert UnitGenerators(DirichletGroup(N)) eq char_gens;
    F := CyclotomicField(order);
    chi := DirichletCharacterFromValuesOnUnitGenerators(DirichletGroup(N,F),[F|F.1^e:e in v]);
    return MinimalBaseRingCharacter(chi);
end function;

function MakeCharacter_3_a_Hecke(Kf)
    return MakeCharacter_3_a();
end function;


function ExtendMultiplicatively(weight, aps, character)
    prec := NextPrime(NthPrime(#aps)) - 1; // we will able to figure out a_0 ... a_prec
    primes := PrimesUpTo(prec);
    prime_powers := primes;
    assert #primes eq #aps;
    log_prec := Floor(Log(prec)/Log(2)); // prec < 2^(log_prec+1)
    F := Universe(aps);
    FXY<X, Y> := PolynomialRing(F, 2);
    // 1/(1 - a_p T + p^(weight - 1) * char(p) T^2) = 1 + a_p T + a_{p^2} T^2 + ...
    R<T> := PowerSeriesRing(FXY : Precision := log_prec + 1);
    recursion := Coefficients(1/(1 - X*T + Y*T^2));
    coeffs := [F!0: i in [1..(prec+1)]];
    coeffs[1] := 1; //a_1
    for i := 1 to #primes do
        p := primes[i];
        coeffs[p] := aps[i];
        b := p^(weight - 1) * F!character(p);
        r := 2;
        p_power := p * p;
        //deals with powers of p
        while p_power le prec do
            Append(~prime_powers, p_power);
            coeffs[p_power] := Evaluate(recursion[r + 1], [aps[i], b]);
            p_power *:= p;
            r +:= 1;
        end while;    
    end for;
    Sort(~prime_powers);
    for pp in prime_powers do
        for k := 1 to Floor(prec/pp) do
            if GCD(k, pp) eq 1 then
                coeffs[pp*k] := coeffs[pp]*coeffs[k];
            end if;
        end for;
    end for;
    return coeffs;
end function;


function qexpCoeffs()
    // To make the coeffs of the qexp of the newform in the Hecke field type "qexpCoeffs();"
    weight := 10;
    raw_aps := [[18], [81], [-1530], [9128], [21132], [31214], [-279342], [144020], [-1763496], [4692510], [-369088], [9347078], [-7226838], [-23147476], [22971888], [78477174], [-20310660], [-179339938], [274528388], [-36342648], [-247089526], [191874800], [-276159276], [-678997350], [-567657502], [1622823462], [-1751025976], [-1542957012], [457665230], [3267943794], [928879328], [-988658748], [-57325302], [-4650524980], [-1402362090], [10154751032], [9366053438], [-734779996], [1555841448], [19044819054], [4328516340], [-956757418], [-17643820608], [13099326914], [4745885958], [-19210229080], [-18127312108], [3910810304], [-35812588572], [5304447110], [-25417154646], [-63349422480], [101928795122], [96474476412], [-7995908862], [106269250824], [-100886359410], [-28634476048], [41225078198], [121800028122], [23282019164], [-91280097306], [-43906990732], [-76483478088], [268488258554], [152055600318], [-263057340148], [297182797778], [140218516188], [-314890729090], [-531194758686], [-497683661400], [277743102608], [255318851414], [81210280700], [-413122928256], [253052022150], [-65356172242], [-240886414638], [-378340063270], [968765698980], [-143788604218], [206889561072], [-48188303566], [437941020440], [758664977724], [-137279666430], [131091223178], [-803507833458], [-1785816198736], [644366062548], [-1877005156320], [1493552300168], [-747313211988], [1012997786420], [203398574904], [593698558590], [-1749501180918], [320554846604], [350228556542], [-2460370242652], [1043558222958], [3818493333684], [484062265530], [17704421372], [-746985322942], [-4339477542132], [1901196887634], [-5496489625320], [4776982081562], [824999503808], [-1469852997946], [5173730633898], [-4802740911700], [5939922195992], [6477383884482], [-1745192195836], [4972442023368], [-1450572111666], [-3051747968940], [-1883147739658], [7843585659554], [2324197335078], [-6892645350996], [11143009010612], [-5780148784578], [-3399020676250], [271271771280], [-2138634975112], [4600332337214], [8602445327780], [1314074790504], [-7958097175408], [-6959363757802], [-13863177783558], [13560448958210], [4967402852934], [-5258091267052], [16181648783838], [-13927645283670], [13818981957932], [3485780849142], [-12147214578856], [-17070225294852], [5580709360670], [12864247581960], [10346474541494], [5632158492378], [-22828883213860], [4629182781024], [5336110351118], [7878899572722], [13897291208564], [-5866069308552], [8522502315788], [16336850963472], [32471231370680], [-12457146318750], [9290118716138], [-32483106483858], [-41625787937292], [18647906982714], [-23298953952472], [-36362793092148], [51814128153618], [-25486233910056], [39283044182912], [49838869596998]];
    aps := ConvertToHeckeField(raw_aps);
    chi := MakeCharacter_3_a_Hecke(Universe(aps));
    return ExtendMultiplicatively(weight, aps, chi);
end function;


// To make the newform (type ModFrm), type "MakeNewformModFrm_3_10_a_b();".
// This may take a long time!  To see verbose output, uncomment the SetVerbose lines below.
// The precision argument determines an initial guess on how many Fourier coefficients to use.
// This guess is increased enough to uniquely determine the newform.
function MakeNewformModFrm_3_10_a_b(:prec:=1)
    chi := MakeCharacter_3_a();
    f_vec := qexpCoeffs();
    Kf := Universe(f_vec);
    // SetVerbose("ModularForms", true);
    // SetVerbose("ModularSymbols", true);
    S := CuspidalSubspace(ModularForms(chi, 10));
    S := BaseChange(S, Kf);
    maxprec := NextPrime(997) - 1;
    while true do
        trunc_vec := Vector(Kf, [0] cat [f_vec[i]: i in [1..prec]]);
        B := Basis(S, prec + 1);
        S_basismat := Matrix([AbsEltseq(g): g in B]);
        if Rank(S_basismat) eq Min(NumberOfRows(S_basismat), NumberOfColumns(S_basismat)) then
            S_basismat := ChangeRing(S_basismat,Kf);
            f_lincom := Solution(S_basismat,trunc_vec);
            f := &+[f_lincom[i]*Basis(S)[i] : i in [1..#Basis(S)]];
            return f;
        end if;
        error if prec eq maxprec, "Unable to distinguish newform within newspace";
        prec := Min(Ceiling(1.25 * prec), maxprec);
    end while;
end function;


// To make the Hecke irreducible modular symbols subspace (type ModSym)
// containing the newform, type "MakeNewformModSym_3_10_a_b();".
// This may take a long time!  To see verbose output, uncomment the SetVerbose line below.
// The default sign is -1.  You can change this with the optional parameter "sign".
function MakeNewformModSym_3_10_a_b( : sign := -1)
    R<x> := PolynomialRing(Rationals());
    chi := MakeCharacter_3_a();
    // SetVerbose("ModularSymbols", true);
    Snew := NewSubspace(CuspidalSubspace(ModularSymbols(chi,10,sign)));
    Vf := Kernel([<2,R![-18, 1]>],Snew);
    return Vf;
end function;