
// Make newform 1620.4.i.h in Magma, downloaded from the LMFDB on 20 April 2024.
// To make the character of type GrpDrchElt, type "MakeCharacter_1620_i();"
// To make the character of type GrpDrchElt with Codomain the HeckeField, type "MakeCharacter_1620_i_Hecke();"
// To make the coeffs of the qexp of the newform in the Hecke field type "qexpCoeffs();"
// To make the newform (type ModFrm), type "MakeNewformModFrm_1620_4_i_h();".
// 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_1620_4_i_h();".
// 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 := [1, -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_1620_i();"
function MakeCharacter_1620_i()
    N := 1620;
    order := 3;
    char_gens := [811, 1541, 1297];
    v := [3, 2, 3];
    // 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;

// To make the character of type GrpDrchElt with Codomain the HeckeField, type "MakeCharacter_1620_i_Hecke();"
function MakeCharacter_1620_i_Hecke(Kf)
    N := 1620;
    order := 3;
    char_gens := [811, 1541, 1297];
    char_values := [[1, 0], [0, -1], [1, 0]];
    assert UnitGenerators(DirichletGroup(N)) eq char_gens;
    values := ConvertToHeckeField(char_values : pass_field := true, Kf := Kf); // the value of chi on the gens as elements in the Hecke field
    F := Universe(values);// the Hecke field
    chi := DirichletCharacterFromValuesOnUnitGenerators(DirichletGroup(N,F),values);
    return chi;
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, 0], [0, 5], [-17, 17], [-30, 30], [0, 61], [120, 0], [-43, 0], [0, 90], [-90, 90], [0, -8], [317, 0], [0, -30], [220, -220], [-180, 180], [630, 0], [0, -840], [-599, 599], [0, -107], [210, 0], [-421, 0], [-353, 353], [-1350, 1350], [-1020, 0], [997, -997], [-960, 960], [0, -1181], [-330, 0], [1454, 0], [0, 1230], [1280, 0], [0, -2220], [-1170, 1170], [0, 1393], [0, 1380], [2659, -2659], [0, -1850], [1121, 0], [0, 300], [-1620, 1620], [630, 0], [-2299, 0], [900, -900], [0, -3461], [-4560, 0], [-2077, 0], [0, 4021], [-80, 80], [-3750, 3750], [0, 1234], [-5880, 0], [0, 5130], [7231, -7231], [-7530, 0], [0, 4560], [-2100, 2100], [-3120, 0], [3449, 0], [-3770, 3770], [-6720, 6720], [0, 100], [0, 870], [3440, 0], [0, -5880], [-1841, 1841], [3420, -3420], [5641, -5641], [0, 1057], [0, -1290], [-3467, 3467], [7740, -7740], [-8130, 0], [-12113, 12113], [0, -4349], [-7663, 0], [0, -3960], [1320, -1320], [14390, 0], [0, 8610], [0, -7097], [0, -8880], [5479, -5479], [12510, 0], [-6790, 0], [11176, -11176], [-13860, 13860], [4740, 0], [1690, -1690], [14700, -14700], [0, 331], [8580, 0], [-14790, 14790], [13097, 0], [0, 1590], [0, -17264], [14730, 0], [0, 870], [6990, 0], [12119, 0], [-21511, 0], [10807, -10807], [20460, 0], [0, -4980], [15840, -15840], [0, 24391], [7673, 0], [18930, -18930], [-14190, 0], [0, -12870], [-19598, 19598], [0, 15163], [-29599, 0], [0, -2490], [-3713, 3713], [19409, 0], [6480, -6480], [0, -30260], [-21510, 0], [0, -540], [29280, -29280], [0, 21769], [14551, -14551], [-1500, 1500], [-16530, 0], [25972, -25972], [10230, 0], [14623, -14623], [-6690, 0], [-11720, 11720], [0, 4750], [-30724, 0], [0, 960], [0, -22781], [32387, 0], [0, 25290], [0, -16283], [31050, 0], [0, 5053], [0, 43680], [-1920, 0], [37268, 0], [5280, -5280], [0, -4511], [10140, 0], [-11923, 0], [21300, -21300], [40771, -40771], [25170, -25170], [0, 42757], [-15300, 0], [0, 23743], [27810, 0], [-42991, 0], [0, -45600], [4993, -4993], [-4920, 4920], [5456, 0], [-14790, 14790], [-14023, 0], [0, -35250], [-37920, 37920], [49080, 0], [0, 33073], [6420, 0], [0, 25110], [52410, -52410], [52619, 0], [53890, -53890], [-3307, 0], [0, 7320], [-6270, 6270], [0, 5038], [0, 18000], [-12491, 12491], [0, 64033], [0, -30150], [-34748, 34748], [44610, 0], [-41641, 0], [-45194, 45194], [17437, -17437], [-45300, 45300], [0, -25490], [18180, 0], [0, -68340], [-32790, 32790], [-20413, 0], [-34169, 34169], [0, 34063], [-4680, 0], [-44269, 0], [-63960, 63960], [60731, 0], [-74460, 74460], [-62400, 0], [0, -22230], [0, -3719], [13129, -13129], [58770, -58770], [-28170, 0], [0, 9900], [-27101, 27101], [0, 69613], [21463, -21463], [76170, 0], [21270, 0], [-60307, 0], [0, 16290], [-10950, 10950], [0, 87259], [46337, 0], [0, -52440], [14869, -14869], [-50370, 50370], [0, 27930], [-23051, 23051], [0, -38660], [-19050, 19050], [14370, 0], [0, 42030], [0, 39391], [0, 37564], [0, 95820], [-58480, 0], [0, -270], [-16463, 16463], [69540, -69540], [31590, 0], [11083, -11083], [-20100, 20100], [0, 13429], [-30796, 0], [0, -100112], [0, 108030], [-68789, 68789], [-13140, 13140], [0, -83693], [-73140, 0], [0, 47550], [-86220, 0], [-92550, 92550], [5381, 0], [0, -92939], [35603, 0], [0, 41340], [-75450, 75450], [31763, 0], [0, -34680], [0, 25516], [94830, 0], [0, 60550], [-9600, 0], [0, -13140], [38287, -38287], [-121350, 121350], [-96900, 0], [63761, 0], [114280, -114280], [-71370, 0], [22970, 0], [-29189, 29189], [-69870, 69870], [0, -117857], [-70291, 0], [0, 46440], [52036, -52036], [-77850, 0], [74310, -74310], [0, 340], [0, 23160], [0, -3851], [31853, 0], [5110, -5110], [0, 77383], [0, -135710], [112721, 0], [0, 88950], [-37826, 37826], [-9298, 0], [-98070, 98070], [0, 30030], [0, 132892], [-141090, 141090], [26569, -26569], [0, -13553], [116160, 0], [-137869, 0], [0, -47760], [-20783, 20783], [58950, 0], [-138150, 138150], [142830, 0], [0, 25170], [0, 79620], [-20990, 20990], [0, 52680], [53701, -53701], [80040, -80040], [89430, 0], [-32483, 32483], [0, -89921], [-12210, 0], [14036, 0]];
    aps := ConvertToHeckeField(raw_aps);
    chi := MakeCharacter_1620_i_Hecke(Universe(aps));
    return ExtendMultiplicatively(weight, aps, chi);
end function;


// To make the newform (type ModFrm), type "MakeNewformModFrm_1620_4_i_h();".
// 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_1620_4_i_h(:prec:=2)
    chi := MakeCharacter_1620_i();
    f_vec := qexpCoeffs();
    Kf := Universe(f_vec);
    // SetVerbose("ModularForms", true);
    // SetVerbose("ModularSymbols", true);
    S := CuspidalSubspace(ModularForms(chi, 4));
    S := BaseChange(S, Kf);
    maxprec := NextPrime(1999) - 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_1620_4_i_h();".
// 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_1620_4_i_h( : sign := -1)
    R<x> := PolynomialRing(Rationals());
    chi := MakeCharacter_1620_i();
    // SetVerbose("ModularSymbols", true);
    Snew := NewSubspace(CuspidalSubspace(ModularSymbols(chi,4,sign)));
    Vf := Kernel([<7,R![289, 17, 1]>,<11,R![900, 30, 1]>],Snew);
    return Vf;
end function;