
// Make newform 1014.4.a.g in Magma, downloaded from the LMFDB on 25 April 2024.
// To make the character of type GrpDrchElt, type "MakeCharacter_1014_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_1014_4_a_g();".
// 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_1014_4_a_g();".
// 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_1014_a();"
function MakeCharacter_1014_a()
    N := 1014;
    order := 1;
    char_gens := [677, 847];
    v := [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_1014_a_Hecke(Kf)
    return MakeCharacter_1014_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 := [[2], [-3], [-6], [16], [-12], [0], [-126], [-20], [168], [30], [88], [-254], [-42], [-52], [96], [198], [660], [-538], [-884], [-792], [-218], [-520], [492], [-810], [-1154], [-618], [128], [-1476], [-1190], [-462], [-2536], [2292], [726], [380], [-1590], [-2432], [614], [1852], [2136], [1758], [-540], [1982], [-2688], [2302], [-4374], [-1600], [3332], [-2648], [-2244], [5650], [4698], [1200], [718], [6012], [-2046], [-6072], [-6930], [-1352], [-1186], [-2442], [2828], [-4758], [8476], [4632], [-4822], [3426], [2788], [434], [6684], [-2630], [7422], [10440], [10424], [3278], [-6140], [3072], [6150], [106], [1758], [3670], [-9660], [-8462], [-9792], [-7342], [10640], [-17412], [1710], [646], [6018], [6712], [5364], [-9840], [-1424], [-4548], [-6500], [12168], [21090], [-5238], [8588], [-3062], [-8476], [12546], [-12], [19290], [-12148], [10366], [-7644], [-8658], [25800], [16202], [-24136], [4642], [6726], [21220], [-29792], [-10158], [-29828], [1944], [26718], [4260], [-22862], [-32542], [14214], [7092], [13228], [28062], [27250], [-14400], [17984], [-16598], [-1460], [30072], [-18088], [24734], [22278], [-16130], [-29718], [-9524], [-33906], [-630], [20788], [43098], [-14272], [-13644], [-2410], [-23160], [-32078], [-14406], [30620], [-17568], [21706], [-14958], [-32812], [-38856], [-28276], [8112], [-26080], [-49170], [48314], [-34782], [25116], [-15462], [736], [-29268], [-16674], [31272], [-15928], [42014], [-12530], [-19722], [18420], [53818], [13992], [28582], [-19240], [-8310], [-39692], [-26982], [-17152], [-31210], [-61544], [24132], [9038], [-28314], [-32208], [-44970], [51734], [47212], [21670], [-29952], [58178], [23172], [55172], [-9102], [1044], [-59142], [5038], [3638], [-1794], [74088], [-14910], [-38968], [38626], [-37150], [-37260], [11934], [-21080], [15108], [59910], [56972], [6194], [-87978], [-17872], [-64524], [-92760], [-61562], [-8296], [49602], [75816], [15762], [-60658], [29680], [-46770], [47032], [59604], [85070], [59178], [-12240], [62624], [95748], [-14362], [10100], [28168], [68442], [-46132], [40416], [37870], [75402], [60], [-12168], [-64308], [-89908], [62768], [44170], [69582], [-13560], [-84584], [-8052], [39940], [-6192], [-63254], [-23202], [69096], [58970], [-55122], [-71460], [542], [-94124], [46134], [119446], [-80872], [62004], [-11630], [86438], [122286], [-15820], [-132030], [-106122], [120292], [27078], [53462], [-14564], [-54218], [-11320], [42194], [-29648], [24996], [33370], [20518], [83172], [82608], [-65792], [-70296], [-13502], [6676], [-97488], [30098], [23466], [125440], [-19710], [-57858], [-107676], [-114378], [70788], [-25162], [66030], [-2248], [45318], [-106020], [40156], [49178], [-81054], [78680], [26628], [-148228], [-110594], [91884], [42950], [-35880], [-161602], [76368], [-120870], [-31758], [-66652], [99576], [-169990], [164220], [-18672], [100222], [-150510], [187652], [-28874], [-147618], [-488], [-21738], [120242], [9340], [-200828], [86304], [13158], [-42662], [9294], [-135880], [43692], [-24572], [-93564], [-148070], [83262], [30598], [-112696], [-99758], [1434], [-145590], [-60928], [-34878], [-78300], [-18898], [-66044], [-62112], [5094], [-200612], [23686], [26418], [26408], [154190], [15498], [-9600], [91908], [-20766], [22632], [-44894], [-194118], [-71856], [47940], [-189196], [217658], [-85134], [-6608], [-41398], [188292], [8300], [-42288], [-123450], [202592], [249014], [184620], [-196128], [-122978], [-143550], [49034], [-17742], [-242538], [-75584], [-245154], [77020], [67032], [190888], [-57346], [61972], [-122736], [-122690], [198042], [139140], [85844], [-199992], [-24982], [248200], [-37590], [-177748], [165738], [-150950], [-56142], [-188264], [72966], [29610], [-47248], [249946], [267198], [-187508], [-164100], [-201742], [-149226], [272988], [-54268], [7994], [44622], [-120960], [-223216], [-139554], [200328], [182190], [-245054], [18096], [61980], [155302], [35346], [-194892], [206550], [83428], [196680]];
    aps := ConvertToHeckeField(raw_aps);
    chi := MakeCharacter_1014_a_Hecke(Universe(aps));
    return ExtendMultiplicatively(weight, aps, chi);
end function;


// To make the newform (type ModFrm), type "MakeNewformModFrm_1014_4_a_g();".
// 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_1014_4_a_g(:prec:=1)
    chi := MakeCharacter_1014_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_1014_4_a_g();".
// 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_1014_4_a_g( : sign := -1)
    R<x> := PolynomialRing(Rationals());
    chi := MakeCharacter_1014_a();
    // SetVerbose("ModularSymbols", true);
    Snew := NewSubspace(CuspidalSubspace(ModularSymbols(chi,4,sign)));
    Vf := Kernel([<5,R![6, 1]>,<7,R![-16, 1]>],Snew);
    return Vf;
end function;