
// Make newform 1386.4.a.j in Magma, downloaded from the LMFDB on 24 April 2024.
// To make the character of type GrpDrchElt, type "MakeCharacter_1386_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_1386_4_a_j();".
// 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_1386_4_a_j();".
// 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_1386_a();"
function MakeCharacter_1386_a()
    N := 1386;
    order := 1;
    char_gens := [155, 199, 1135];
    v := [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_1386_a_Hecke(Kf)
    return MakeCharacter_1386_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], [0], [-2], [-7], [-11], [26], [46], [-48], [128], [146], [-128], [-26], [-10], [52], [544], [-318], [48], [466], [516], [392], [754], [0], [-624], [1590], [1018], [-474], [-984], [-92], [1246], [1630], [1016], [1920], [-1482], [-2608], [-1310], [-192], [658], [2236], [1664], [662], [2540], [2762], [16], [5138], [-4350], [-4040], [-1820], [-2360], [6416], [-1558], [-522], [-2152], [-606], [1608], [4446], [-6600], [1854], [-272], [-5010], [-314], [3480], [4230], [-1552], [4864], [4786], [1530], [-10844], [402], [5980], [-3094], [4494], [2752], [-2024], [5246], [3892], [-6752], [-12486], [1938], [-4530], [-13718], [-15280], [478], [-6280], [13802], [8728], [3540], [-4194], [-14134], [-234], [13696], [-16104], [11272], [-304], [-10572], [-15004], [-16872], [-818], [2270], [12776], [-23050], [-6564], [4522], [20440], [16518], [-8828], [-15550], [4536], [11142], [16248], [-9646], [4064], [-15098], [4470], [-21184], [-8760], [3582], [-23168], [-30216], [-8158], [-11932], [26882], [13090], [33790], [-24588], [-1128], [18786], [12102], [-18112], [12728], [17138], [-8340], [-8304], [2152], [-25594], [17190], [-15086], [-14178], [-18304], [38206], [-7146], [-21256], [-38670], [-21112], [3172], [30346], [9480], [-24958], [26806], [23128], [-12496], [30478], [-25506], [-13244], [-25456], [-51652], [46392], [17832], [41334], [23058], [11678], [36436], [-21098], [7184], [40048], [-51938], [28968], [37504], [20298], [-51614], [55750], [44052], [18098], [-4408], [18154], [39632], [-54250], [-21508], [-30710], [-54072], [-4430], [-52632], [-44600], [17398], [15662], [-6328], [-11754], [-5930], [-28248], [45242], [61160], [-55454], [-34348], [43220], [-10162], [-7012], [30390], [56978], [73982], [63078], [-5256], [25274], [16648], [70546], [68186], [59904], [84054], [5120], [62148], [47382], [28616], [48674], [9566], [9104], [25216], [69816], [-59534], [-86368], [-54498], [1176], [26770], [89334], [98392], [95166], [-90584], [-88032], [77254], [38534], [-41976], [73912], [8716], [15422], [-48424], [93760], [-2298], [-6288], [10800], [112410], [-1662], [-54876], [-95120], [-13932], [93640], [-37440], [17030], [-82842], [-98784], [-76904], [-36328], [56052], [-96096], [78006], [-78218], [67832], [-31430], [98878], [9804], [55478], [2704], [-96962], [43898], [30088], [-91996], [-80694], [-55582], [38934], [30576], [-4614], [56966], [76156], [-19622], [-70670], [-39476], [-104886], [86072], [19594], [-44752], [-55620], [-134746], [-78998], [-129840], [109632], [-95416], [92880], [92554], [-53816], [10200], [134770], [-134254], [-76824], [-22842], [-119430], [-2744], [-22506], [-150248], [-165194], [47830], [-80376], [-120202], [-6608], [-165528], [-119014], [-15934], [22976], [-81276], [-134188], [68306], [-118180], [-161174], [121536], [125598], [76848], [-40790], [-14786], [-61956], [-148104], [97554], [-113200], [-29568], [81426], [182718], [-16840], [39578], [24782], [-112160], [30294], [20634], [181292], [-116832], [-101400], [-38070], [105750], [-173198], [-120760], [-39560], [25996], [-124944], [16894], [-50394], [134442], [-141768], [-42994], [6006], [188734], [-146408], [81426], [-52524], [-27814], [40708], [-41712], [19662], [-72024], [19706], [-148294], [-133288], [136206], [50622], [-16336], [-88888], [-6802], [106128], [23918], [145790], [47000], [18324], [12680], [-158838], [153686], [-114560], [-184390], [-98820], [-245656], [107072], [-160190], [94904], [41238], [140712], [186040], [-105238], [-20066], [-156678], [-89570], [102326], [171224], [213694], [37336], [-19168], [-216304], [104658], [276460], [157600], [134578], [-262642], [-21964], [199472], [83472], [32682], [270488], [-110930], [-118036], [268906], [106866], [-200034], [75120], [-135474], [-129138], [-31808], [-169610], [31710], [-45376], [79184], [259338], [-151974], [-36612], [-188492], [100522], [21534], [117856], [-11608], [-114010], [192808], [260450], [240954], [-190384], [-154384], [150338], [126054], [222612], [234630], [-56968], [111968]];
    aps := ConvertToHeckeField(raw_aps);
    chi := MakeCharacter_1386_a_Hecke(Universe(aps));
    return ExtendMultiplicatively(weight, aps, chi);
end function;


// To make the newform (type ModFrm), type "MakeNewformModFrm_1386_4_a_j();".
// 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_1386_4_a_j(:prec:=1)
    chi := MakeCharacter_1386_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_1386_4_a_j();".
// 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_1386_4_a_j( : sign := -1)
    R<x> := PolynomialRing(Rationals());
    chi := MakeCharacter_1386_a();
    // SetVerbose("ModularSymbols", true);
    Snew := NewSubspace(CuspidalSubspace(ModularSymbols(chi,4,sign)));
    Vf := Kernel([<5,R![2, 1]>,<13,R![-26, 1]>],Snew);
    return Vf;
end function;