
// Make newform 1800.4.a.t in Magma, downloaded from the LMFDB on 18 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_t();".
// 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_t();".
// 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], [2], [-39], [84], [61], [151], [58], [-192], [-18], [-138], [-229], [-164], [212], [-578], [336], [858], [-209], [780], [-403], [-230], [1293], [1369], [382], [794], [-1348], [775], [446], [231], [2386], [-2452], [1125], [-1377], [-1920], [1854], [-634], [103], [-44], [1128], [2245], [3050], [4222], [-3357], [166], [3372], [5601], [828], [2388], [-2844], [-5962], [4320], [3857], [287], [-2130], [3066], [3744], [-3346], [7040], [3010], [-6001], [-4802], [6149], [878], [-4042], [-3844], [-2717], [-1603], [11607], [4030], [-2106], [-7394], [-6940], [7486], [1285], [-9622], [-1974], [-8084], [5667], [-4835], [-4619], [7476], [-7810], [-2029], [3208], [-13227], [-3617], [6215], [7108], [-3364], [18964], [10926], [-4350], [1324], [9068], [19836], [-2682], [3035], [7701], [-18112], [19541], [-13508], [8712], [9623], [604], [3629], [-9219], [19111], [17086], [9035], [-14784], [-17846], [-11618], [-9556], [-19394], [-12138], [27036], [17556], [17262], [-10517], [1408], [-9626], [-28464], [3963], [-31781], [28004], [-35228], [8658], [5728], [21460], [29164], [-29478], [576], [2880], [-20789], [26421], [-32504], [996], [-15134], [36942], [-11748], [1198], [-6788], [-33011], [17732], [-8480], [30014], [-21643], [2799], [19384], [5132], [-4430], [24317], [26100], [24356], [29900], [-34838], [-26334], [30949], [-25276], [1216], [6033], [-41792], [2105], [-30119], [-18438], [-2230], [-6804], [-21499], [-38904], [-17821], [21338], [8206], [-9907], [-35224], [-9351], [48009], [-46278], [-14472], [28874], [-13756], [-12707], [50090], [-21561], [-52068], [-25460], [-37466], [9217], [8209], [-44112], [-37519], [45657], [16513], [64074], [24431], [-47398], [-22205], [-274], [56091], [-68486], [-35258], [77852], [50624], [35555], [3652], [-5214], [-11074], [-35108], [24882], [20560], [-39842], [-50450], [20576], [76715], [-87518], [-21191], [56018], [-75229], [-88238], [-35718], [65078], [-10804], [60159], [-74738], [36883], [3598], [-32049], [16416], [-51388], [92355], [-98546], [79709], [-90802], [-36583], [-63984], [-109182], [53125], [76872], [-62665], [-43828], [21747], [62717], [-64998], [-93820], [-31805], [-11904], [45974], [29484], [15765], [9656], [118784], [120357], [-71054], [47435], [-107672], [22251], [1600], [40811], [-13062], [97843], [64430], [-14551], [4024], [-64652], [94925], [-3821], [9466], [-60094], [55877], [63218], [18644], [18603], [-82399], [96328], [-38370], [-99300], [47307], [146904], [-30930], [24021], [-116676], [84922], [-65768], [57460], [25127], [78120], [25026], [33968], [-90416], [148062], [-130210], [-75372], [24627], [-40561], [-145900], [-326], [41326], [-3672], [101340], [124291], [10703], [3446], [-134192], [-33753], [-66327], [60799], [-34924], [13470], [-70366], [-14486], [-41398], [114198], [167215], [-57616], [-22636], [-74330], [-165991], [-36504], [-35343], [-42051], [-122596], [-37413], [-58620], [-96534], [51610], [-414], [35329], [69340], [104262], [109556], [47764], [103394], [68890], [158711], [-95741], [-167776], [-10632], [-186750], [-73673], [-178254], [-134200], [158937], [-17686], [-42844], [-66194], [-54244], [162644], [51623], [-89776], [-135162], [100445], [-179986], [-109590], [30748], [111194], [-70325], [198468], [12493], [-130418], [-22700], [-77930], [885], [-168510], [-199076], [-206437], [-74579], [-48888], [-207080], [-217267], [146028], [-60001], [-3212], [-74572], [58592], [-16654], [-164655], [-156390], [50613], [87951], [-66098], [64638], [-72647], [-51288], [-14505], [-187624], [54494], [-50752], [9696], [-214936], [28646], [76683], [270580], [170471], [83123], [27050], [16443], [83116], [150247], [-185617], [149606], [-54440], [-197139], [221724], [45905], [236640], [-143144], [-176314], [-49482], [216668], [102951], [8193], [-202608], [77461], [-26660], [199969], [235442], [301614], [152962], [-79815], [-199302], [-278720], [-127606], [-125826], [48643], [-224557], [-43344], [-308240], [-62979], [-1613], [-130914]];
    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_t();".
// 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_t(: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_t();".
// 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_t( : sign := -1)
    R<x> := PolynomialRing(Rationals());
    chi := MakeCharacter_1800_a();
    // SetVerbose("ModularSymbols", true);
    Snew := NewSubspace(CuspidalSubspace(ModularSymbols(chi,4,sign)));
    Vf := Kernel([<7,R![-2, 1]>,<11,R![39, 1]>,<17,R![-61, 1]>],Snew);
    return Vf;
end function;