/* Copyright (C) 2021 Alessandro Languasco */ /**************** A. LANGUASCO ******************** ************* COMPUTATION OF Gamma_T, quadratic case *******/ \\ Global variables: global(tab); \\global variable used to initialise sumnum global(defaultprecision); {gammaT_2331(Pbound,defaultprecision)=local(minutes, millisec, seconds, coeff, quadform1, quadform2, elaptimefinalcomp, s1, s2, s3, w, qminusone, S, psquared, pcubed, ptoqminusone, ptoq, isquadpmodq, logp, errs1, errs2, errs3, errS, Pbound1, derlog, D, gammaK, Resultfile, gammaT, s2quad, s2quadbis, gammaquad, gammaquadbis ); \\ minutes,millisec,seconds: used just to compute the elapsed computation time; local variables \\ defaultprecision: used to fix the precision used in the computations; global variable print("************ A. LANGUASCO *************"); print("********* COMPUTATION of gammaT_quad **********"); print("******* q = 23 or q = 31 ********"); \\ precision setting default(realprecision,defaultprecision+10); default(format,f); Pbound1=Pbound+1; tab = sumnuminit(); \\ needed for the sumnum functions in the errors Resultfile = fileopen("Resultfile2331.csv", "w"); filewrite(Resultfile,"q;r=(q-1)/2;gammaKr(q);gamma_quad(q);err_gammaquad;gammaquad-err;gammaquad+err;gammaT_quad(q);err_gammaT;gammaT-err;gammaT+err"); gettime(); print("Starting computation of gammaT"); print("internal sums over primes p up to ", Pbound); print("---------"); forprimestep(q=23,31, Mod(3,4), errs1=0; errs2=0; errs3=0; s1=0; s2=0; s3=0; s2quad = 0; s2quadbis = 0; S=0; qminusone = q-1; w = (q+1)/2; quadform1 = Qfb(2,1,w/4); \\ create the quadratic form 2*x^2 + xy + w/4*y^2 quadform2 = Qfb(1,1,w/2); \\ create the quadratic form 1*x^2 + xy + w/2*y^2 forprime(p=2, qminusone, \\ to keep q=3 in this case logp = log(p); psquared = p*p*1.0; pcubed = psquared*p; ptoqminusone = p^qminusone*1.0; ptoq = ptoqminusone * p; isquadpmodq = kronecker(p,q); \\ p^((q-1)/2) mod q if( isquadpmodq == 1, s2quadbis += ( q/(ptoq -1) - qminusone/(ptoqminusone - 1) ) * logp ; if (qfbsolve(quadform1,p) <> [], s2 += ( 2/(psquared - 1) - 3/(pcubed - 1) ) * logp; s2quad += ( 2/(psquared - 1) - 3/(pcubed - 1) + q/(ptoq -1) - qminusone/(ptoqminusone - 1) ) * logp ); \\ p is in S2 of the paper if (qfbsolve(quadform2,p) <> [], s3 += ( q/(ptoq -1) - qminusone/(ptoqminusone - 1) ) * logp) \\ p is in S3 of the paper , s1 += logp/(psquared - 1); ); ); forprime(p=qminusone+3, Pbound, \\ from q+2 logp = log(p); psquared = p*p*1.0; pcubed = psquared*p; ptoqminusone = p^qminusone*1.0; ptoq = ptoqminusone * p; isquadpmodq = kronecker(p,q); \\ p^((q-1)/2) mod q if( isquadpmodq == 1, s2quadbis += ( q/(ptoq -1) - qminusone/(ptoqminusone - 1) ) * logp ; if (qfbsolve(quadform1,p) <> [], s2 += ( 2/(psquared - 1) - 3/(pcubed - 1) ) * logp; s2quad += ( 2/(psquared - 1) - 3/(pcubed - 1) + q/(ptoq -1) - qminusone/(ptoqminusone - 1) ) * logp ); \\ p is in S2 of the paper if (qfbsolve(quadform2,p) <> [], s3 += ( q/(ptoq -1) - qminusone/(ptoqminusone - 1) ) * logp) \\ p is in S3 of the paper , s1 += logp/(psquared - 1); ); ); S = s2 - s1 - s3; \\ REMARK : just errs3 depends on q and becomes tiny VERY quickly errs1 = sumnum(n=Pbound1, log(n)*( 2/(n^2 -1) - 3/(n^3 -1) ), tab); \\ at most one time I have the contribution of primes s.t. the summand is log(n)*( 2/(n^2 -1) - 3/(n^3 -1) errs2 = sumnum(n=Pbound1, log(n)/(n^2 -1) , tab); errs3 = sumnum(n=Pbound1, log(n)*( q/(n^q -1) - qminusone/(n^qminusone -1) ), tab); \\print(errs1); \\print(errs2); \\print(errs3); errS = abs(errs1) + abs(errs2) + abs(errs3); D = kronecker(-1,q)*q; \\q = 3 mod 4 \\L = lfuncreate(D); \\ create the quadratic L-function mod q \\derlog = lfun(L,1,1)/lfun(L,1); \\ log derivative at 1 for L derlog = Lquadquot(1, D); \\ log derivative at 1 for the quadratic L-function mod q gammaK = Euler + derlog; coeff = 0.5 * ( gammaK - log(q)/qminusone); \\ coefficient in the sum gammaT = S + coeff ; gammaquad = gammaT - s2quad; \\ using displayed formula of the paper before remark 7 gammaquadbis = coeff - s1 - s2quadbis ; \\ using formula (42) of the paper \\print("Sr = ", - s1 - s2quadbis); print("q = ",q); print("gammaK(",qminusone/2,",",q,") = ", gammaK ); print("gammaTquad(",q,") = ", gammaT ); print("VIA GAMMAT --- gammaquad(",qminusone/2,",",q,") = ", gammaquad ); print("DIRECTLY ----- gammaquad(",qminusone/2,",",q,") = ", gammaquadbis ); print("err(",q,") = ", errS); print("***** err affects ONLY gamma_quad and gammaT_quad; NOT gammaK"); print("***** its q-depending part becomes TINY very QUICKLY"); if (abs(S) < 10*errS, print("********* recompute !!!")); print("---------"); filewrite1(Resultfile,q); filewrite1(Resultfile,";"); filewrite1(Resultfile,qminusone/2); filewrite1(Resultfile,";"); filewrite1(Resultfile, gammaK); filewrite1(Resultfile,";"); filewrite1(Resultfile, gammaquad); filewrite1(Resultfile,";"); filewrite1(Resultfile, errS); filewrite1(Resultfile,";"); filewrite1(Resultfile, gammaquad-errS); filewrite1(Resultfile,";"); filewrite1(Resultfile, gammaquad+errS); filewrite1(Resultfile,";"); filewrite1(Resultfile, gammaT); filewrite1(Resultfile,";"); filewrite1(Resultfile, errS); filewrite1(Resultfile,";"); filewrite1(Resultfile, gammaT-errS); filewrite1(Resultfile,";"); filewrite(Resultfile, gammaT+errS); fileflush(Resultfile); ); fileclose(Resultfile); elaptimefinalcomp=gettime(); seconds=floor(elaptimefinalcomp/1000)%60; minutes=floor(elaptimefinalcomp/60000); millisec=elaptimefinalcomp- minutes*60000 - seconds*1000; print("gammaT(q) computation time (output time included): ", minutes, " min, ", seconds, " sec, ", millisec, " millisec"); print("****** END PROGRAM ********"); } /* quadratic L function attached to the Kronecker symbol for small |D| */ {Lquad(s, D)= local(chi, q); print(D); chi=vector(abs(D),i,kronecker(D,i)); q=length(chi); print(q); if (s==1, -sum(a=1,q-1,chi[a]*psi(a/q))/q \\ L(1,chi) , sum(a=1,q-1,chi[a]*zetahurwitz(s,a/q))/q^s \\ L(s,chi) ) }; /* log derivative quadratic L function attached to the Kronecker symbol for small |D| */ {Lquadquot(s, D)= local(chi, q); chi=vector(abs(D),i,kronecker(D,i)); q=length(chi); if (s==1, Euler + log(2*Pi) + q*sum(a=1,q-1, chi[a]* lngamma(a/q)) / sum(a=1,q-1, chi[a]* a) \\ L'/L(1,chi) , -log(q)+ sum(a=1,q-1,chi[a]* zetahurwitz(s,a/q,1)) /sum(a=1,q-1,chi[a]* zetahurwitz(s,a/q)) \\ L'/L(s,chi) ) }; /**** RESULTS **** gp2.13.1 and gp2c0.0.12 compiled by myself, see below languasc@languasc-macmini modular_forms % gp2c-run -pmy_ -g -W new-gammaT_2331.gp Warning:new-gammaT_2331.gp:35: variable undeclared f Reading GPRC: /Users/languasc/.gprc GPRC Done. GP/PARI CALCULATOR Version 2.13.2 (released) i386 running darwin (x86-64/GMP-6.2.1 kernel) 64-bit version compiled: Jun 24 2021, Apple clang version 12.0.5 (clang-1205.0.22.9) threading engine: single (readline v8.0 enabled, extended help enabled) Copyright (C) 2000-2020 The PARI Group PARI/GP is free software, covered by the GNU General Public License, and comes WITHOUT ANY WARRANTY WHATSOEVER. Type ? for help, \q to quit. Type ?17 for how to get moral (and possibly technical) support. parisizemax = 2048000000, primelimit = 500000 ------------------ 10^7 macmini ------------------ ? gammaT_2331(10^7,19) ************ A. LANGUASCO ************* ********* COMPUTATION of gammaT_quad ********** ******* q = 23 or q = 31 ******** Starting computation of gammaT internal sums over primes p up to 10000000 --------- q = 23 gammaK(11,23) = 0.15509471686308553509479869689 gammaTquad(23) = 0.21669143132172819553707663173 VIA GAMMAT --- gammaquad(11,23) = -0.15489093208491103537727577947 DIRECTLY ----- gammaquad(11,23) = -0.15489093208491103537727577947 err(23) = 0.0000051354282042446748453048730480 ***** err affects ONLY gamma_quad and gammaT_quad; NOT gammaK ***** its q-depending part becomes TINY very QUICKLY --------- q = 31 gammaK(15,31) = 0.12605851277979231152685968549 gammaTquad(31) = 0.15610590636741710245343806162 VIA GAMMAT --- gammaquad(15,31) = -0.20075888743896822864578785578 DIRECTLY ----- gammaquad(15,31) = -0.20075888743896822864578785578 err(31) = 0.0000051354282042446748453048730480 ***** err affects ONLY gamma_quad and gammaT_quad; NOT gammaK ***** its q-depending part becomes TINY very QUICKLY --------- gammaT(q) computation time (output time included): 0 min, 10 sec, 489 millisec ****** END PROGRAM ******** ------------------ 10^8 macmini ------------------ ? gammaT_2331(10^8,19) ************ A. LANGUASCO ************* ********* COMPUTATION of gammaT_quad ********** ******* q = 23 or q = 31 ******** Starting computation of gammaT internal sums over primes p up to 100000000 --------- q = 23 gammaK(11,23) = 0.15509471686308553509479869689 gammaTquad(23) = 0.21669144632953378943498233416 VIA GAMMAT --- gammaquad(11,23) = -0.15489097707906845381843813714 DIRECTLY ----- gammaquad(11,23) = -0.15489097707906845381843813714 err(23) = 0.00000058262041671736679632398346071 ***** err affects ONLY gamma_quad and gammaT_quad; NOT gammaK ***** its q-depending part becomes TINY very QUICKLY --------- q = 31 gammaK(15,31) = 0.12605851277979231152685968549 gammaTquad(31) = 0.15610592139737023942382226963 VIA GAMMAT --- gammaquad(15,31) = -0.20075893243723704390789582759 DIRECTLY ----- gammaquad(15,31) = -0.20075893243723704390789582759 err(31) = 0.00000058262041671736679632398346071 ***** err affects ONLY gamma_quad and gammaT_quad; NOT gammaK ***** its q-depending part becomes TINY very QUICKLY --------- gammaT(q) computation time (output time included): 1 min, 34 sec, 559 millisec ****** END PROGRAM ******** ------------------ 10^9 macmini ------------------ ? gammaT_2331(10^9,19) ************ A. LANGUASCO ************* ********* COMPUTATION of gammaT_quad ********** ******* q = 23 or q = 31 ******** Starting computation of gammaT internal sums over primes p up to 1000000000 --------- q = 23 gammaK(11,23) = 0.15509471686308553509479869689 gammaTquad(23) = 0.21669144783001254375317383925 VIA GAMMAT --- gammaquad(11,23) = -0.15489098157892463761322734809 DIRECTLY ----- gammaquad(11,23) = -0.15489098157892463761322734809 err(23) = 0.000000065169797447919436019899666606 ***** err affects ONLY gamma_quad and gammaT_quad; NOT gammaK ***** its q-depending part becomes TINY very QUICKLY --------- q = 31 gammaK(15,31) = 0.12605851277979231152685968549 gammaTquad(31) = 0.15610592289513110477765660215 VIA GAMMAT --- gammaquad(15,31) = -0.20075893693785241938125323158 DIRECTLY ----- gammaquad(15,31) = -0.20075893693785241938125323158 err(31) = 0.000000065169797447919436019899666606 ***** err affects ONLY gamma_quad and gammaT_quad; NOT gammaK ***** its q-depending part becomes TINY very QUICKLY --------- gammaT(q) computation time (output time included): 14 min, 15 sec, 79 millisec ****** END PROGRAM ******** ------------------ 10^10 macmini ------------------ ? gammaT_2331(10^(10),19) ************ A. LANGUASCO ************* ********* COMPUTATION of gammaT_quad ********** ******* q = 23 or q = 31 ******** Starting computation of gammaT internal sums over primes p up to 10000000000 --------- q = 23 gammaK(11,23) = 0.15509471686308553509479869689 gammaTquad(23) = 0.21669144798002462200016437079 VIA GAMMAT --- gammaquad(11,23) = -0.15489098202892158539149576929 DIRECTLY ----- gammaquad(11,23) = -0.15489098202892158539149576929 err(23) = 0.0000000072077552782838615242249215450 ***** err affects ONLY gamma_quad and gammaT_quad; NOT gammaK ***** its q-depending part becomes TINY very QUICKLY --------- q = 31 gammaK(15,31) = 0.12605851277979231152685968549 gammaTquad(31) = 0.15610592304511585773529332629 VIA GAMMAT --- gammaquad(15,31) = -0.20075893738785799911664718597 DIRECTLY ----- gammaquad(15,31) = -0.20075893738785799911664718597 err(31) = 0.0000000072077552782838615242249215450 ***** err affects ONLY gamma_quad and gammaT_quad; NOT gammaK ***** its q-depending part becomes TINY very QUICKLY --------- gammaT(q) computation time (output time included): 138 min, 39 sec, 31 millisec ****** END PROGRAM ******** *******/