
// Make newform 720.4.a.r in Magma, downloaded from the LMFDB on 10 May 2024.
// To make the character of type GrpDrchElt, type "MakeCharacter_720_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_720_4_a_r();".
// 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_720_4_a_r();".
// 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_720_a();"
function MakeCharacter_720_a()
    N := 720;
    order := 1;
    char_gens := [271, 181, 641, 577];
    v := [1, 1, 1, 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_720_a_Hecke(Kf)
    return MakeCharacter_720_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 := 4;
    raw_aps := [[0], [0], [5], [-20], [-24], [74], [-54], [124], [-120], [78], [-200], [-70], [-330], [-92], [-24], [-450], [24], [-322], [196], [-288], [-430], [520], [156], [-1026], [-286], [1734], [-452], [-1404], [-1474], [-1086], [-1244], [2328], [-2118], [-2324], [-258], [808], [2378], [52], [-3720], [-426], [-1440], [-3130], [3576], [2666], [2718], [3832], [-1100], [-1964], [660], [-1906], [1458], [1176], [866], [432], [-2526], [5448], [2574], [3184], [3962], [8286], [2716], [-6018], [-9236], [1536], [-7342], [3894], [-3692], [-8998], [5244], [6302], [-3414], [4824], [3508], [10802], [-1460], [-4872], [14046], [-2734], [15942], [8714], [11976], [11054], [720], [-15622], [9880], [-16116], [-9018], [-3670], [-17562], [-1172], [6876], [2280], [3076], [-18912], [-9956], [-10656], [2766], [-10530], [-12692], [18110], [-3620], [14166], [-13404], [18654], [7684], [-1726], [10596], [-2862], [-23592], [-9574], [-17444], [-2374], [12162], [-8804], [12688], [9150], [-25292], [-2736], [-22218], [14520], [-10618], [1370], [13758], [11988], [-32996], [25902], [-27394], [34848], [-28028], [18002], [-15284], [-18768], [-8696], [-38662], [-23874], [23618], [-11538], [14884], [11334], [-44730], [42748], [31686], [-11036], [25884], [15950], [13800], [-27862], [7314], [28780], [-32688], [36650], [2646], [-10892], [-43464], [14884], [-1248], [6640], [-29946], [45002], [-6090], [56388], [-10854], [42316], [24480], [6906], [6960], [-47792], [9938], [-2590], [7374], [-23064], [-22282], [-7176], [40898], [47608], [41214], [-38900], [5454], [-9452], [30710], [-44660], [-5784], [-18934], [32634], [29112], [33630], [-15766], [46276], [63722], [-47016], [-29878], [62796], [-4988], [-5466], [17244], [-31686], [-77614], [26426], [-37326], [35568], [-16458], [13912], [44642], [-14254], [-1512], [-18810], [9304], [-75612], [-498], [71476], [-40822], [-59442], [-18740], [-29076], [-84576], [32618], [-75548], [-80058], [-52848], [-7962], [-29266], [-36632], [-83490], [31804], [87756], [6590], [85722], [-58464], [-57836], [-36456], [84170], [-34076], [67984], [22782], [26044], [-16656], [55442], [-41778], [52296], [-26352], [-99564], [-18524], [40372], [88238], [-63798], [28704], [-15212], [-93360], [-115004], [57888], [-83950], [-85122], [-41880], [-102166], [61062], [-12168], [-32938], [-28244], [56718], [110162], [-2252], [-67020], [-31138], [-125854], [36066], [126412], [85734], [38022], [-125300], [106014], [72614], [23884], [46658], [-21512], [-61414], [40012], [144492], [102350], [506], [-31560], [-23520], [-28424], [12720], [-108754], [-17516], [-109080], [79058], [99822], [63448], [18270], [42054], [6876], [-25278], [56904], [-31102], [57270], [141928], [-39954], [-49320], [117052], [49826], [-29898], [26368], [36636], [168532], [144434], [101316], [-11914], [13056], [-74878], [-8448], [-103626], [-34410], [135052], [-136944], [152426], [10272], [123432], [122330], [102558], [-85700], [149330], [414], [-24212], [30498], [-27934], [-86780], [8692], [-2472], [22014], [-117178], [15366], [128680], [-59316], [-3188], [177588], [-48922], [18714], [52058], [162916], [1442], [-25230], [-27834], [-86312], [31254], [-161280], [-21922], [-144812], [-89880], [125838], [-41996], [-93238], [-2850], [193156], [-27034], [-86526], [92976], [181944], [-117582], [-75624], [-114766], [-116634], [220392], [-233736], [243844], [-33190], [-71994], [-55316], [21386], [83424], [-204356], [50280], [-80346], [-87464], [152210], [-92064], [80928], [125714], [-24138], [-82294], [73878], [-58182], [145708], [36186], [-41132], [-255528], [-103304], [-67174], [-152732], [215376], [-171646], [-128226], [11592], [-28460], [200568], [20378], [-88400], [-96474], [-56540], [130926], [-147562], [175986], [-241580], [166890], [90774], [-60320], [-24046], [-222906], [5812], [-272616], [290594], [-176466], [-186300], [-58916], [172946], [-26514], [-271584], [235564], [-98910], [134664], [-7338], [98642], [-80808], [126384], [-258766], [-255786], [11100], [36630], [121444], [152184]];
    aps := ConvertToHeckeField(raw_aps);
    chi := MakeCharacter_720_a_Hecke(Universe(aps));
    return ExtendMultiplicatively(weight, aps, chi);
end function;


// To make the newform (type ModFrm), type "MakeNewformModFrm_720_4_a_r();".
// 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_720_4_a_r(:prec:=1)
    chi := MakeCharacter_720_a();
    f_vec := qexpCoeffs();
    Kf := Universe(f_vec);
    // SetVerbose("ModularForms", true);
    // SetVerbose("ModularSymbols", true);
    S := CuspidalSubspace(ModularForms(chi, 4));
    S := BaseChange(S, Kf);
    maxprec := NextPrime(2999) - 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_720_4_a_r();".
// 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_720_4_a_r( : sign := -1)
    R<x> := PolynomialRing(Rationals());
    chi := MakeCharacter_720_a();
    // SetVerbose("ModularSymbols", true);
    Snew := NewSubspace(CuspidalSubspace(ModularSymbols(chi,4,sign)));
    Vf := Kernel([<7,R![20, 1]>,<11,R![24, 1]>],Snew);
    return Vf;
end function;