Skip to content

Commit

Permalink
added trace on subspace of newforms. Works for T_7*W on S_2(77).
Browse files Browse the repository at this point in the history
  • Loading branch information
@assaferan
assaferan committed Jun 22, 2023
1 parent 6545c97 commit df9f600
Show file tree
Hide file tree
Showing 4 changed files with 157 additions and 32 deletions.
52 changes: 39 additions & 13 deletions src/traceAL.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -3,24 +3,50 @@

#include "traceAL_tools.h"

int main()
int main(int argc, char* argv[])
{
GEN /* f, NK,*/ k, N, al, upto, from, prec;
// long prec; // , rem;
GEN al;
long k, N, upto, from, prec;
int only_primes, newspace;

if ((argc < 5) || (argc > 8))
{
printf("Incorrect number of arguments.\n");
printf("Usage: %s <level> <weight> <num_traces> <new> (<only_primes>) (<from> <upto>) \n", argv[0]);
printf("computes traces of S_k(N)^+ for level N between from and to, weight k and for Hecke operators up to the greatest between the Sturm bound, 30 sqrt(p) and 1000 if not specified.\n");
return -1;
}
N = atoi(argv[1]);
k = atoi(argv[2]);
prec = atoi(argv[3]);
newspace = atoi(argv[4]);
newspace = (newspace ? 1 : 0); // making sure it is either 1 or 0
if (newspace)
only_primes = 1;
else
only_primes = atoi(argv[5]);
if (argc >= 6 + (1-newspace)) {
from = atoi(argv[5 + (1 - newspace)]);
upto = atoi(argv[6 + (1 - newspace)]);
}

pari_init(10000000000,2);
printf("N = "); N = gp_read_stream(stdin);
printf("k = "); k = gp_read_stream(stdin);
printf("prec = "); prec = gp_read_stream(stdin);

time_t start = time(NULL);
al = traceALupto(gtos(N), gtos(k), gtos(prec));
if (only_primes)
al = traceALprimes(N, k, prec, newspace);
else
al = traceALupto(N, k, prec);

printf("single run took %ld seconds\n", time(NULL)-start);
pari_printf("al = %Ps\n", al);
// NK = mkvec2(N,k);
// f = mftraceform(NK,0);
printf("from = "); from = gp_read_stream(stdin);
printf("upto = "); upto = gp_read_stream(stdin);
time_t timing = timeTraceAL(gtos(upto), gtos(from), gtos(k), -1, 1);
printf("took %ld seconds\n", timing);

if (argc >= 6 + (1-newspace)) {
from = atoi(argv[5 + (1 - newspace)]);
upto = atoi(argv[6 + (1 - newspace)]);
time_t timing = timeTraceAL(upto, from, k, -1, only_primes, newspace);
printf("took %ld seconds\n", timing);
}
pari_close_mf();
pari_close();
return 0;
Expand Down
113 changes: 102 additions & 11 deletions src/traceAL_tools.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -540,37 +540,128 @@ GEN traceAL(long N, long n, long k)
return ret;
}

// At the moment only workd for primes
GEN traceALNewTrivial(long N, long k)
{
GEN trace = gen_0;
GEN NN = mkintn(1,N);
GEN div_N = divisors(NN);
long num_divs_N = lg(div_N);

for (long idx = 1; idx < num_divs_N; idx++) {
GEN D = gel(div_N, idx);
GEN D2 = gmul(D,D);

// if ((N % (d*d)) != 0) continue;
if (gmod(NN, D2) != gen_0) continue;
long N_div_d2 = gtos(gdivent(NN, D2));
trace = gadd(trace, gmulsg(moebius(D),traceAL(N_div_d2, 1, k)));
}

// pari_printf("Trace of W_{%Ps} on new subspace is: %Ps\n", NN, trace);

return trace;
}

GEN traceALNewTrivialContribution(long N, long p, long k)
{
GEN trace;
GEN trace2 = gen_0;
GEN trace3 = gen_0;
GEN div_N = divisors(mkintn(1,N));
long num_divs_N = lg(div_N);

for (long idx = 1; idx < num_divs_N; idx++) {
long d = gtos(gel(div_N, idx));
long N_div_d = N / d;
// We use valuation, since p^3 might not be a single-word integer
if ((d % p == 0) || (gvaluation(mkintn(1,N_div_d), mkintn(1,p)) >= 3)) continue;
if (! issquare(mkintn(1,p*N_div_d)) ) continue;
trace2 = gadd(trace2, traceALNewTrivial(d, k));
}
// pari_printf("Contribution from N2 for level %d for prime %d and weight %d is: %Ps \n", N, p, k, trace2);

if (N % p == 0) {
long N_div_p = N / p;
GEN div_N_div_p = divisors(mkintn(1,N_div_p));
long num_divs_N_div_p = lg(div_N_div_p);
for (long idx = 1; idx < num_divs_N_div_p; idx++) {
long d = gtos(gel(div_N_div_p, idx));
if (! issquare(mkintn(1,N_div_p / d)) ) continue;
trace3 = gadd(trace3, traceALNewTrivial(d, k));
}
}
// pari_printf("Contribution from N3 for level %d for prime %d and weight %d is: %Ps \n", N, p, k, trace3);

// Since p is at most one word length, we can estimate the size of its powers accordingly
trace = gsub(gmul(gpow(mkintn(1,p), mkintn(1,k/2), k/2), trace3),
gmul(gpow(mkintn(1,p), mkintn(1,k/2-1), k/2-1), trace2));

// pari_printf("Trivial contribution from level %d for prime %d and weight %d is: %Ps \n", N, p, k, trace);
return trace;
}

// At the moment only works for primes
GEN traceALNew(long N, long p, long k)
{
GEN ret = gen_0;
GEN trace = gen_0;
GEN NN = mkintn(1,N);
GEN div_N = divisors(NN);
long num_divs_N = lg(div_N);

return ret;
// printf("In traceALNew with N = %ld.\n", N);

for (long idx = 1; idx < num_divs_N; idx++) {
GEN D = gel(div_N, idx);
GEN D2 = gmul(D,D);

// pari_printf("D = %Ps \n", D);
// pari_printf("N mod (D^2) = %Ps \n", gmod(NN,D2));
// printf("(N mod (D^2) != 0) = %d \n", (gmod(NN,D2) != 0));
// printf("gen_0 (N mod (D^2) != 0) = %d \n", (gmod(NN,D2) != gen_0));
// printf("cmp (N mod (D^2) != 0) = %d \n", gcmp(gmod(NN,D2),gen_0));
// printf("gtos(D) mod p = %ld\n", gtos(D) % p);
// printf("(gtos(D) mod p == 0) = %d\n", gtos(D) % p == 0);

if (gmod(NN, D2) != gen_0) continue;
if (gtos(D) % p == 0) continue;

long N_div_d2 = gtos(gdivent(NN, D2));
trace = gadd(trace, gmulsg(moebius(D), gsub(traceAL(N_div_d2,p,k), traceALNewTrivialContribution(N_div_d2,p,k))));
// pari_printf("Accumulated trace is: %Ps \n", trace);
}

return trace;
}

GEN traceALprimes(long N, long k, long prec)
GEN traceALprimes(long N, long k, long prec, int newspace)
{
GEN p_list = primes0(mkvec2(mkintn(1,1), nextprime(mkintn(1,prec))));
long num_primes = lg(p_list);
GEN res = cgetg(num_primes-1, t_VEC);
long p;

GEN (*trace_func)(long, long, long);
trace_func = (newspace ? &traceALNew : &traceAL);

// printf("In traceALprimes. num_primes = %ld. newspace = %d. \n", num_primes, newspace);
// printf("trace_func = %x.\n", (unsigned int)trace_func);

for (long idx = 1; idx < num_primes - 1; idx++)
{
p = gtos(gel(p_list, idx));
gel(res, idx) = traceAL(N, p, k);
gel(res, idx) = (*trace_func)(N, p, k);
}
return res;
}

GEN trace_primes(long N, long k, long prec)
GEN trace_primes(long N, long k, long prec, int newspace)
{
GEN p_list = primes0(mkvec2(mkintn(1,1), nextprime(mkintn(1,prec))));
long num_primes = lg(p_list);
GEN res = cgetg(num_primes-1, t_VEC);
long p;
GEN NK = mkvec2(mkintn(1,N),mkintn(1,k));
GEN f = mftraceform(NK,1);
GEN f = mftraceform(NK,newspace ? 0 : 1);

for (long idx = 1; idx < num_primes - 1; idx++)
{
Expand All @@ -591,7 +682,7 @@ GEN traceALupto(long N, long k, long prec)
return res;
}

time_t timeTraceAL(long upTo, long from, long k, long num_traces, int only_primes)
time_t timeTraceAL(long upTo, long from, long k, long num_traces, int only_primes, int newspace)
{
time_t start = time(NULL);
long prec;
Expand Down Expand Up @@ -620,12 +711,12 @@ time_t timeTraceAL(long upTo, long from, long k, long num_traces, int only_prime

// printf("p = %ld, prec = %ld\n", p, prec);
if (only_primes)
res = traceALprimes(N, k, prec+1);
res = traceALprimes(N, k, prec+1, newspace);
else
res = traceALupto(N, k, prec+1);

if (only_primes)
coefs = trace_primes(N, k, prec+1);
coefs = trace_primes(N, k, prec+1, newspace);
else {
NK = mkvec2(mkintn(1,N),mkintn(1,k));
f = mftraceform(NK,1);
Expand All @@ -641,6 +732,6 @@ time_t timeTraceAL(long upTo, long from, long k, long num_traces, int only_prime
//pari_printf("traces := %Ps;\ntracesAL := %Ps;\n", coefs, res);
fclose(outfile);
}
printf("Finished.\n");
// printf("Finished.\n");
return time(NULL) - start;
}
6 changes: 3 additions & 3 deletions src/traceAL_tools.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -35,10 +35,10 @@ long alpha(ulong n);

GEN traceAL(long N, long n, long k);

GEN traceALprimes(long N, long k, long prec);
GEN traceALprimes(long N, long k, long prec, int newspace);

GEN trace_primes(long N, long k, long prec);
GEN trace_primes(long N, long k, long prec, int newspace);

GEN traceALupto(long N, long k, long prec);

time_t timeTraceAL(long upTo, long from, long k, long num_traces, int only_primes);
time_t timeTraceAL(long upTo, long from, long k, long num_traces, int only_primes, int newspace);
18 changes: 13 additions & 5 deletions src/traceALbatch.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -6,10 +6,10 @@
int
main(int argc, char* argv[])
{
if ((argc < 4) || (argc > 6))
if ((argc < 4) || (argc > 7))
{
printf("Incorrect number of arguments.\n");
printf("Usage: %s <from> <to> <weight> (<num_traces>) (<only_primes>) \n", argv[0]);
printf("Usage: %s <from> <to> <weight> (<num_traces>) (<only_primes>) (<new>) \n", argv[0]);
printf("computes traces of S_k(N)^+ for level N between from and to, weight k and for Hecke operators up to the greatest between the Sturm bound, 30 sqrt(p) and 1000 if not specified.\n");
return -1;
}
Expand All @@ -18,14 +18,22 @@ main(int argc, char* argv[])
long k = atoi(argv[3]);
long num_traces = -1;
int only_primes = 0;
int newspace = 0;

if (argc >= 5)
num_traces = atoi(argv[4]);
if (argc == 6)
if (argc >= 6)
only_primes = atoi(argv[5]);
if (argc == 7)
newspace = atoi(argv[6]);

if ((newspace) && (!only_primes)) {
printf("Currently only supports newspace traces for Hecke operators at primes.");
return -1;
}

pari_init(10000000000,2);
timeTraceAL(upto, from, k, num_traces, only_primes);
// time_t timing = timeTraceAL(upto, from, k);
timeTraceAL(upto, from, k, num_traces, only_primes, newspace);
// printf("took %ld seconds\n", timing);
pari_close_mf();
pari_close();
Expand Down

0 comments on commit df9f600

Please sign in to comment.