
// Make newform 648.4.i.v in Magma, downloaded from the LMFDB on 19 April 2024.
// To make the character of type GrpDrchElt, type "MakeCharacter_648_i();"
// To make the character of type GrpDrchElt with Codomain the HeckeField, type "MakeCharacter_648_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_648_4_i_v();".
// 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_648_4_i_v();".
// 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 := [169, 0, -104, 0, 51, 0, -8, 0, 1];
        Kf := NumberField(Polynomial([elt : elt in poly]));
        AssignNames(~Kf, ["nu"]);
    end if;
    Rf_num := [[1, 0, 0, 0, 0, 0, 0, 0], [-169, 0, 408, 0, -51, 0, 8, 0], [1300, 5980, 0, -1632, 0, 204, 13, -32], [1300, -9542, 0, 2754, 0, -510, 13, 28], [100, 0, 0, 0, 0, 0, 1, 0], [-676, 1976, 969, 1632, -204, -204, 32, 32], [1976, 0, -969, 0, 204, 0, -19, 0], [-676, -2236, 969, -2754, -204, 510, 32, -106]];
    Rf_basisdens := [1, 663, 663, 663, 17, 663, 221, 663];
    Rf_basisnums := ChangeUniverse([[z : z in elt] : elt in Rf_num], Kf);
    Rfbasis := [Rf_basisnums[i]/Rf_basisdens[i] : i in [1..Degree(Kf)]];
    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_648_i();"
function MakeCharacter_648_i()
    N := 648;
    order := 3;
    char_gens := [487, 325, 569];
    v := [3, 3, 2];
    // 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_648_i_Hecke();"
function MakeCharacter_648_i_Hecke(Kf)
    N := 648;
    order := 3;
    char_gens := [487, 325, 569];
    char_values := [[1, 0, 0, 0, 0, 0, 0, 0], [1, 0, 0, 0, 0, 0, 0, 0], [0, -1, 0, 0, 0, 0, 0, 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, 0, 0, 0], [0, 0, 0, 0, 0, 0, 0, 0], [0, 2, 0, 0, 0, 0, 0, 1], [0, 0, -1, -1, 0, 1, 0, 1], [2, -2, 2, 1, 0, -2, -1, -1], [0, -1, 0, 0, 2, -1, -2, -3], [-4, 0, -4, 2, 1, 0, 0, 0], [20, 0, -5, 3, 0, 0, 0, 0], [0, 50, 0, 0, 7, -6, -7, 5], [54, -54, -4, -9, 0, 4, 4, 9], [0, -20, 0, 0, 16, 2, -16, 2], [-69, 0, 17, -1, -22, 0, 0, 0], [0, 96, 0, 0, -20, 4, 20, -12], [-40, 40, 3, 19, 0, -3, -56, -19], [192, -192, 0, 16, 0, 0, -12, -16], [-236, 0, 24, 10, -34, 0, 0, 0], [0, 248, 0, 0, 28, 24, -28, -12], [137, -137, 1, 3, 0, -1, -114, -3], [0, -116, 0, 0, 96, -21, -96, 19], [-430, 0, -6, -27, 43, 0, 0, 0], [-191, 0, -40, -44, -36, 0, 0, 0], [-172, 172, -23, -23, 0, 23, -160, 23], [532, -532, -20, -38, 0, 20, -6, 38], [-528, 0, -56, -40, -55, 0, 0, 0], [454, -454, 4, 40, 0, -4, -204, -40], [420, -420, 60, 54, 0, -60, 26, -54], [0, -556, 0, 0, 112, 46, -112, -98], [-276, 0, 4, 94, 134, 0, 0, 0], [-287, 0, 7, 121, -74, 0, 0, 0], [0, 912, 0, 0, -171, -40, 171, 72], [172, 0, 103, -57, -128, 0, 0, 0], [0, 418, 0, 0, 243, -86, -243, -87], [828, -828, -84, -70, 0, 84, 125, 70], [0, -828, 0, 0, 120, -42, -120, 110], [0, -110, 0, 0, -342, 52, 342, 13], [-440, 440, 66, 178, 0, -66, -128, -178], [0, 797, 0, 0, 34, -29, -34, -79], [1440, 0, -96, 48, 32, 0, 0, 0], [0, 450, 0, 0, 419, 90, -419, -3], [-102, 102, -64, 21, 0, 64, 174, -21], [-192, 0, -8, -72, 308, 0, 0, 0], [-1278, 0, -36, -244, 80, 0, 0, 0], [438, -438, 46, 247, 0, -46, -367, -247], [0, 1549, 0, 0, -156, 58, 156, 166], [-6, 0, -44, -95, -296, 0, 0, 0], [748, 0, -112, -128, 128, 0, 0, 0], [0, -1576, 0, 0, -88, -183, 88, -71], [-268, 268, -199, -87, 0, 199, 560, 87], [-1482, 1482, -34, -125, 0, 34, 201, 125], [0, 851, 0, 0, -334, -47, 334, -9], [1192, 0, 68, 336, -327, 0, 0, 0], [0, -1148, 0, 0, 698, -92, -698, 186], [207, -207, 290, -6, 0, -290, 416, 6], [1382, 0, -58, 187, 437, 0, 0, 0], [0, -2696, 0, 0, -269, -196, 269, -40], [-2776, 2776, -248, -124, 0, 248, -504, 124], [4082, 0, -172, 45, -50, 0, 0, 0], [-772, 0, 119, 263, 256, 0, 0, 0], [278, -278, 240, -72, 0, -240, 680, 72], [-2984, 2984, 148, 144, 0, -148, -15, -144], [0, -2036, 0, 0, -584, 196, 584, 124], [0, -1722, 0, 0, -166, -188, 166, 55], [2652, 0, 70, -146, 184, 0, 0, 0], [0, -4428, 0, 0, 606, 4, -606, 98], [1413, -1413, -458, -82, 0, 458, 912, 82], [-3258, 3258, -216, -29, 0, 216, 178, 29], [-2612, 2612, 133, 189, 0, -133, 464, -189], [0, 126, 0, 0, -388, -420, 388, 40], [0, -1710, 0, 0, -901, 554, 901, 209], [-42, 42, -232, 364, 0, 232, 708, -364], [2608, -2608, 72, 112, 0, -72, 48, -112], [-918, 0, 290, 529, -645, 0, 0, 0], [2672, -2672, 206, 30, 0, -206, 304, -30], [0, -590, 0, 0, -192, 0, 192, 576], [-172, 0, 34, -742, -152, 0, 0, 0], [0, 2898, 0, 0, 411, -246, -411, -987], [1124, -1124, -408, 514, 0, 408, 6, -514], [1619, 0, -133, 229, -178, 0, 0, 0], [0, -1644, 0, 0, 385, 464, -385, 370], [0, -345, 0, 0, 28, -6, -28, -130], [0, -1716, 0, 0, -1390, -452, 1390, -26], [-2537, 2537, 611, -11, 0, -611, -922, 11], [6452, 0, -428, 258, -814, 0, 0, 0], [2167, 0, 186, -342, -152, 0, 0, 0], [4644, -4644, 272, 304, 0, -272, -800, -304], [616, -616, -368, -204, 0, 368, 1912, 204], [-4856, 0, -260, 104, 1112, 0, 0, 0], [-5921, 5921, 38, 362, 0, -38, -900, -362], [7100, -7100, 144, -546, 0, -144, -214, 546], [0, 4216, 0, 0, 688, -256, -688, -32], [-6684, 0, -300, 266, 230, 0, 0, 0], [-1970, 1970, -42, -493, 0, 42, -1467, 493], [-9932, 0, -462, 226, -448, 0, 0, 0], [0, 9570, 0, 0, -2297, -206, 2297, 609], [0, -960, 0, 0, 664, -147, -664, -755], [-7692, 0, 100, 450, -1578, 0, 0, 0], [0, 3636, 0, 0, 1798, -392, -1798, -478], [-6472, 0, -60, 504, 552, 0, 0, 0], [6688, 0, -221, 291, -600, 0, 0, 0], [4965, 0, -267, -981, -142, 0, 0, 0], [-1460, 1460, -444, 524, 0, 444, -1352, -524], [-9650, 0, -96, 759, 1238, 0, 0, 0], [0, 5048, 0, 0, -3000, -64, 3000, 444], [13604, -13604, -248, -334, 0, 248, -561, 334], [0, -6816, 0, 0, 2112, -108, -2112, -892], [-1933, 0, -292, -1072, -972, 0, 0, 0], [938, -938, 234, 477, 0, -234, 3123, -477], [6880, 0, -400, -536, 1861, 0, 0, 0], [0, 12228, 0, 0, -374, 964, 374, -494], [-1111, 1111, 100, 1456, 0, -100, -1620, -1456], [0, 16552, 0, 0, 128, -241, -128, 271], [5678, 0, 316, 28, -448, 0, 0, 0], [0, 3024, 0, 0, 1379, -620, -1379, -1100], [3728, -3728, 494, 1390, 0, -494, -1232, -1390], [-7160, 0, 262, -298, 64, 0, 0, 0], [-5236, 5236, -1436, 506, 0, 1436, 849, -506], [0, -1572, 0, 0, 552, 146, -552, 730], [6182, 0, 862, 487, -1863, 0, 0, 0], [0, 11844, 0, 0, 1150, 592, -1150, -214], [-1746, 1746, 150, -297, 0, -150, 3037, 297], [0, -7979, 0, 0, 2206, -725, -2206, 453], [-8011, 8011, 954, 206, 0, -954, -1060, -206], [18524, -18524, -256, 1182, 0, 256, -1286, -1182], [-8312, 0, 336, 1300, -40, 0, 0, 0], [-12236, 12236, -367, -1031, 0, 367, -640, 1031], [17802, 0, 856, 717, 1358, 0, 0, 0], [-2237, 2237, 5, -837, 0, -5, 82, 837], [-4082, 0, 54, 1523, -1895, 0, 0, 0], [18432, -18432, 537, 537, 0, -537, 1552, -537], [0, 8062, 0, 0, -664, -684, 664, -1604], [16748, 0, -48, 24, 328, 0, 0, 0], [0, -13888, 0, 0, 1780, -416, -1780, 1048], [0, 9216, 0, 0, 896, 29, -896, -1699], [13970, 0, -380, 708, -608, 0, 0, 0], [0, 3700, 0, 0, -7, 212, 7, 742], [0, -563, 0, 0, -1356, 34, 1356, 1070], [-5838, 0, -608, -1023, -544, 0, 0, 0], [0, -13744, 0, 0, -3080, 491, 3080, 1163], [0, -8998, 0, 0, 778, -828, -778, 809], [22148, 0, -872, 706, -489, 0, 0, 0], [5512, 0, -644, 92, 3776, 0, 0, 0], [-17682, 17682, -808, 351, 0, 808, -158, -351], [0, 13444, 0, 0, -1968, 80, 1968, -2288], [-2734, 0, 146, 545, 2407, 0, 0, 0], [-4622, 0, -916, -924, 1208, 0, 0, 0], [14526, -14526, 246, -373, 0, -246, -4207, 373], [4810, -4810, 1648, 824, 0, -1648, 3112, -824], [-17892, 17892, 1256, -26, 0, -1256, -1671, 26], [0, -9564, 0, 0, -2184, -1628, 2184, -116], [24358, 0, 206, -105, -859, 0, 0, 0], [0, -3003, 0, 0, -4966, 861, 4966, 827], [224, 0, 1756, -492, 44, 0, 0, 0], [10436, 0, -378, 254, 376, 0, 0, 0], [0, -8822, 0, 0, 3147, 1154, -3147, -327], [-18160, 18160, -238, 850, 0, 238, 4304, -850], [-28150, 28150, -414, -375, 0, 414, -3781, 375], [-13076, 0, 1079, 2119, -1120, 0, 0, 0], [-11924, 11924, 140, 778, 0, -140, 2107, -778], [-3321, 0, 1788, -2184, 316, 0, 0, 0], [0, -20790, 0, 0, -2058, -104, 2058, -867], [14490, -14490, -1262, -2115, 0, 1262, 3551, 2115], [-2292, 0, 1080, 1062, -3845, 0, 0, 0], [0, 3592, 0, 0, -5840, 1755, 5840, -853], [37062, 0, 702, 475, 2609, 0, 0, 0], [0, 15720, 0, 0, -2984, 668, 2984, -2536], [-5148, 5148, 612, 2154, 0, -612, -4922, -2154], [-1112, 0, -765, -1469, 4688, 0, 0, 0], [-21103, 21103, 681, -1617, 0, -681, 2658, 1617]];
    aps := ConvertToHeckeField(raw_aps);
    chi := MakeCharacter_648_i_Hecke(Universe(aps));
    return ExtendMultiplicatively(weight, aps, chi);
end function;


// To make the newform (type ModFrm), type "MakeNewformModFrm_648_4_i_v();".
// 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_648_4_i_v(:prec:=8)
    chi := MakeCharacter_648_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(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_648_4_i_v();".
// 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_648_4_i_v( : sign := -1)
    R<x> := PolynomialRing(Rationals());
    chi := MakeCharacter_648_i();
    // SetVerbose("ModularSymbols", true);
    Snew := NewSubspace(CuspidalSubspace(ModularSymbols(chi,4,sign)));
    Vf := Kernel([<5,R![48455521, 2617336, 1686718, 27904, 45331, 1024, 286, -8, 1]>],Snew);
    return Vf;
end function;