
// Make newform 230.4.a.f in Magma, downloaded from the LMFDB on 14 May 2024.
// To make the character of type GrpDrchElt, type "MakeCharacter_230_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_230_4_a_f();".
// 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_230_4_a_f();".
// 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
        poly := [-18, -1, 1];
        Kf := NumberField(Polynomial([elt : elt in poly]));
        AssignNames(~Kf, ["nu"]);
    end if;
    Rfbasis := [Kf.1^i : i in [0..Degree(Kf)-1]];
    inp_vec := Vector(Rfbasis)*ChangeRing(Transpose(Matrix([[elt : elt in row] : row in input])),Kf);
    return Eltseq(inp_vec);
end function;


// To make the character of type GrpDrchElt, type "MakeCharacter_230_a();"
function MakeCharacter_230_a()
    N := 230;
    order := 1;
    char_gens := [47, 51];
    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_230_a_Hecke(Kf)
    return MakeCharacter_230_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], [-1, -1], [5, 0], [-8, -1], [4, 10], [-9, 9], [-34, -11], [12, 10], [-23, 0], [-55, -2], [-33, -26], [-242, -7], [-129, -74], [-166, -22], [-297, -5], [-156, 7], [126, 105], [-92, 12], [-286, 41], [679, -104], [-183, 23], [-16, -56], [578, 117], [562, -18], [-1038, -164], [600, -87], [220, -308], [-66, 185], [-202, -12], [972, -17], [-1647, -111], [261, -391], [-1496, -50], [161, -242], [838, -390], [-627, -143], [-1768, 25], [1793, 455], [240, 680], [-366, 546], [7, -505], [-324, 972], [942, -602], [709, 753], [-507, -1033], [88, 206], [-962, -291], [-964, -310], [1446, -94], [-2510, -330], [223, -273], [4051, -78], [-1142, -416], [1020, -546], [-145, 1507], [4174, 171], [-159, 1048], [-514, -433], [773, 979], [542, -632], [3104, -1039], [-1678, -613], [3472, -502], [-1751, 2073], [-8290, 23], [-4594, 286], [1529, -1036], [2786, -1036], [580, -524], [-4539, 1598], [2393, -575], [-9648, 652], [466, 711], [986, 1288], [-5120, -1140], [3660, -695], [-5060, 2354], [-2077, 123], [3774, -324], [-1195, 152], [-6578, -1756], [-5238, -1172], [2166, -82], [2064, 563], [-397, 1819], [595, 2735], [-6796, 1025], [-2130, -1661], [-2677, 2493], [96, 586], [-516, 3623], [-4816, 724], [-4779, 3269], [13749, 1038], [-2641, -248], [-9728, -271], [-15677, -199], [10808, -1886], [-11810, 2682], [-17137, 129], [1089, 2729], [3656, 19], [-94, 2375], [-6834, 1824], [11484, 2540], [15461, 233], [-19417, -1695], [-306, 4674], [4200, -3736], [-6079, -1666], [11536, -3104], [11060, -1182], [16282, -2821], [10020, 2648], [7202, -2418], [19450, 1608], [14780, 1271], [-273, -579], [-8789, 2637], [10524, 2466], [5502, 5248], [6615, -1285], [21944, 2413], [3033, -3283], [1864, -3548], [-11850, -978], [6442, -2022], [-22326, -871], [-12520, 3309], [-6786, -3855], [-7839, -2424], [-3740, -1636], [-12266, -3740], [6682, -6889], [9341, -3016], [41286, -30], [-14874, 3764], [5966, 5281], [6446, 4411], [-4680, 609], [-14897, 3118], [-15574, 24], [-18757, -2301], [-14962, -3865], [40768, 1055], [-14978, 4256], [6806, -2640], [-16017, -1359], [17739, -3698], [-31143, -573], [-1454, 6892], [-26444, 2228], [-25020, 5408], [-48273, -955], [-19242, -2587], [-3372, 1904], [-24302, -5926], [18603, 6816], [-27480, 4642], [12956, 6068], [-22837, -1957], [16298, -8124], [21193, -1095], [41798, 3478], [18810, 2879], [-4932, -3093], [39838, -1783], [12694, -554]];
    aps := ConvertToHeckeField(raw_aps);
    chi := MakeCharacter_230_a_Hecke(Universe(aps));
    return ExtendMultiplicatively(weight, aps, chi);
end function;


// To make the newform (type ModFrm), type "MakeNewformModFrm_230_4_a_f();".
// 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_230_4_a_f(:prec:=2)
    chi := MakeCharacter_230_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(997) - 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_230_4_a_f();".
// 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_230_4_a_f( : sign := -1)
    R<x> := PolynomialRing(Rationals());
    chi := MakeCharacter_230_a();
    // SetVerbose("ModularSymbols", true);
    Snew := NewSubspace(CuspidalSubspace(ModularSymbols(chi,4,sign)));
    Vf := Kernel([<3,R![-16, 3, 1]>],Snew);
    return Vf;
end function;