
// Make newform 1800.4.a.bf in Magma, downloaded from the LMFDB on 25 April 2024.
// To make the character of type GrpDrchElt, type "MakeCharacter_1800_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_1800_4_a_bf();".
// 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_1800_4_a_bf();".
// 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_1800_a();"
function MakeCharacter_1800_a()
    N := 1800;
    order := 1;
    char_gens := [1351, 901, 1001, 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_1800_a_Hecke(Kf)
    return MakeCharacter_1800_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], [0], [19], [-22], [-1], [-58], [-53], [58], [-22], [-35], [270], [468], [431], [-230], [0], [-446], [127], [811], [-36], [-522], [1368], [-1138], [-144], [1079], [1440], [-124], [-432], [701], [-1044], [-504], [-72], [-1850], [1836], [2398], [1871], [3293], [-883], [4104], [1706], [-662], [4121], [958], [-3187], [-2282], [1043], [4139], [-413], [5652], [-4391], [-2052], [4320], [4265], [2412], [6948], [-3564], [1534], [7704], [-4877], [3758], [-935], [6214], [3905], [-9662], [-5147], [-9216], [2196], [-2521], [7034], [7362], [-9382], [10116], [3277], [-10891], [-2591], [-612], [-3708], [3833], [9288], [-755], [9576], [9414], [2254], [6301], [-3779], [7632], [-2988], [1370], [7992], [3096], [-8338], [-4918], [19637], [-3096], [6875], [11268], [-16078], [17366], [4913], [17605], [14560], [-20614], [3442], [-22082], [451], [-20987], [3888], [-11268], [3996], [-24965], [-4176], [9558], [-9000], [-15625], [-31175], [6732], [6228], [-396], [-4702], [22140], [-13518], [11250], [-15170], [22680], [128], [-25682], [-4951], [22190], [-7685], [-29574], [-32580], [-3060], [-7992], [-22841], [17172], [-30869], [34884], [-5039], [540], [41328], [-12853], [-29470], [24407], [15048], [-28406], [26914], [-21275], [39132], [-448], [-26856], [3653], [6552], [-4481], [-14666], [-51804], [-31198], [-27001], [22694], [-29503], [-14566], [42322], [43416], [-21528], [-27050], [-24934], [8388], [29033], [25326], [125], [-13846], [-46670], [29107], [-58910], [-11573], [6443], [-8280], [307], [38110], [1763], [-27413], [2071], [47592], [-18233], [-13622], [-62302], [-4550], [34903], [15949], [-28171], [49082], [3655], [30204], [125], [-35546], [-48074], [-34546], [-45954], [-13806], [14594], [-35302], [-10404], [-10820], [39166], [73890], [11268], [-18202], [76051], [-16034], [50076], [23093], [58877], [16870], [-53153], [-46620], [-30874], [68382], [-9109], [-13910], [-92556], [8482], [51461], [-12113], [-102514], [-52128], [38124], [62299], [32126], [106006], [3709], [55210], [50867], [-59436], [95399], [-94774], [81468], [67176], [-27846], [70344], [-104026], [74016], [-48650], [108685], [-98947], [39511], [-63850], [77890], [-89839], [92146], [-46115], [-29002], [-23166], [-27526], [-37836], [39581], [-19714], [112572], [47789], [-5004], [66600], [56754], [7055], [-62950], [69678], [18325], [82634], [-82208], [134654], [-59040], [-62460], [-85838], [62963], [11072], [-141066], [-4159], [-30329], [-41993], [122566], [-61669], [100871], [121176], [-144828], [-21600], [-91642], [-133433], [-72487], [-45130], [-22625], [146916], [-47304], [72684], [112550], [-97812], [-31932], [-81072], [-20287], [-114804], [-67627], [-612], [-103918], [25164], [165114], [-157068], [142721], [118346], [-28332], [5651], [48852], [160741], [68184], [93634], [-115178], [143676], [-56556], [77075], [-166558], [-114263], [-104638], [-88524], [-130033], [-163210], [141176], [36469], [183528], [38123], [-139334], [9378], [70723], [146197], [-113746], [24012], [-159211], [56542], [15697], [123538], [-8837], [8820], [-179425], [189936], [-59291], [67464], [199153], [-100908], [-15646], [-101341], [-59054], [82130], [-186587], [15821], [156384], [-42444], [-28925], [193878], [88358], [174313], [-193429], [197374], [-176076], [172922], [-61466], [-130990], [-112355], [-166990], [42494], [-80280], [-140868], [-83827], [-181080], [-192188], [-23977], [-214236], [181637], [-111298], [-134820], [-171792], [27631], [-73620], [-257234], [7705], [182556], [-162845], [26626], [-100238], [-49248], [-60660], [-114660], [179410], [69013], [-16110], [146484], [140710], [15517], [-116748], [133632], [-62533], [-18986], [173863], [260477], [-82296], [186173], [-166356], [-172998], [-95904], [-226243], [-155664], [143208], [-162864], [-25201], [138708], [-224227], [237132], [-267751], [257004], [158738], [-6821], [-235333], [5364], [-5746], [94805], [-85068], [-49464], [-32832], [-160399], [177862], [-275530], [57762], [1332], [219780], [223488], [-40339], [137074]];
    aps := ConvertToHeckeField(raw_aps);
    chi := MakeCharacter_1800_a_Hecke(Universe(aps));
    return ExtendMultiplicatively(weight, aps, chi);
end function;


// To make the newform (type ModFrm), type "MakeNewformModFrm_1800_4_a_bf();".
// 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_1800_4_a_bf(:prec:=1)
    chi := MakeCharacter_1800_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_1800_4_a_bf();".
// 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_1800_4_a_bf( : sign := -1)
    R<x> := PolynomialRing(Rationals());
    chi := MakeCharacter_1800_a();
    // SetVerbose("ModularSymbols", true);
    Snew := NewSubspace(CuspidalSubspace(ModularSymbols(chi,4,sign)));
    Vf := Kernel([<7,R![-19, 1]>,<11,R![22, 1]>,<17,R![58, 1]>],Snew);
    return Vf;
end function;