/* Copyright (C) 2021 Alessandro Languasco */ /**************** A. LANGUASCO ******************** *******/ /* Here we compute the values of loggamma(a/q) */ /***************** Precomputations for the external C program ****************/ \\ \p needed to decide the precision default (format, f) for the floating point {ale_loggamma(x,y)= local(vec, aoverq, S, oneminusaoverq, minusaoverq, q, g, a, u , coeff, r1, m, maxindex, datavector, minutes, millisec, abslogaoverq, abslogoneminusaoverq, seconds, elaptimeprecomp, aleloggammafile, valueatonehalf, name, namefile, start, correction); default(format , f); \\ print numbers with floating point notation (not E notation) vec = readvec("./primroot.res"); \\getting q and g from primroot.res q=vec[1]; g=vec[2]; print("Precomputation loggamma-values for q = ",q, " and saving on file"); if (x<0, x=0); \\ needed to avoid problem with the first interval if (x>=y, error("first exponent should be less than the second")); if (y>q-2, y=q-2); \\ needed to avoid problem with the last interval name="./precomp_aleloggamma"; namefile = Strprintf("%s-%012d-%012d.res", name, x, y); aleloggammafile = fileopen(namefile, "w"); if (x==0,filewrite(aleloggammafile,q)); \\ print(q);print(g); m=128; \\ binary precision required for our algorithm coeff = (m+1)*log(2)*1.0; \\ decimal precision \\print(maxindex); valueatonehalf = 0.5*log(Pi); maxindex = m+10; datavector = vector(maxindex,i, zeta(i+1)/(i+1)) ; \\initialisation of the needed coefficients gettime(); a=g^x%q; for(k=x, y, aoverq = (a/q)*1.0; \\ already sorted for the final summation in external S=0; \\ \\if (aoverq == 0.5, print("1/2"); S = valueatonehalf, \\ q is prime => a/q <> 1/2 \\ a/q <>1/2 minusaoverq = -aoverq; oneminusaoverq = 1 - aoverq; abslogaoverq = abs(log(aoverq)); abslogoneminusaoverq = abs(log(oneminusaoverq)); if (aoverq > 0.5, \\ a/q >1/2 r1 = ceil( ( coeff + abslogaoverq ) / abslogoneminusaoverq ); \\ number of summands start = oneminusaoverq; correction =Euler*oneminusaoverq , \\ a/q <1/2 r1 = ceil( ( coeff + abslogoneminusaoverq ) / abslogaoverq ); \\ number of summands start = minusaoverq; correction = abslogaoverq + Euler * minusaoverq; ); u = sqr(start); \\ starting point for computing for (k=2, r1, S+= u * datavector[k-1]; u*= start; ); S = correction+ S; \\ correction values filewrite(aleloggammafile, S); a=(a*g)%q; \\a=g^k%q; ); elaptimeprecomp=gettime(); fileclose(aleloggammafile); /* print computation time */ seconds=floor(elaptimeprecomp/1000)%60; minutes=floor(elaptimeprecomp/60000); millisec=elaptimeprecomp- minutes*60000 - seconds*1000; \\print(elaptimeprecomp); print("Precomputation time (I/O included): ", minutes, " min, ", seconds, " sec, ", millisec, " millisec"); } /************************************ 2.13.0 su optiplex gp2c-run -pmy_ -g -W ale_loggamma.gp ? ale_loggamma(0,10007) Precomputation loggamma-values for q = 10007 and saving on file Precomputation time (I/O included): 0 min, 0 sec, 151 millisec ? preclngamma(0,10007) Precomputation lngamma-values for q = 10007 and saving on file Precomputation time (I/O included): 0 min, 0 sec, 252 millisec --------- ? ale_loggamma(0,305741) Precomputation loggamma-values for q = 305741 and saving on file Precomputation time (I/O included): 0 min, 4 sec, 532 millisec ? preclngamma(0,305741) Precomputation lngamma-values for q = 305741 and saving on file Precomputation time (I/O included): 0 min, 7 sec, 539 millisec --------- ? ale_loggamma(0,6766811) Precomputation loggamma-values for q = 6766811 and saving on file Precomputation time (I/O included): 1 min, 42 sec, 599 millisec ? preclngamma(0,6766811) Precomputation lngamma-values for q = 6766811 and saving on file Precomputation time (I/O included): 2 min, 49 sec, 627 millisec --------- ? ale_loggamma(0,10000019) Precomputation loggamma-values for q = 10000019 and saving on file Precomputation time (I/O included): 2 min, 33 sec, 980 millisec ? preclngamma(0,10000019) Precomputation lngamma-values for q = 10000019 and saving on file Precomputation time (I/O included): 4 min, 10 sec, 632 millisec --------- ? ale_loggamma(0,28227761) Precomputation loggamma-values for q = 28227761 and saving on file Precomputation time (I/O included): 7 min, 16 sec, 158 millisec ? preclngamma(0,28227761) Precomputation lngamma-values for q = 28227761 and saving on file Precomputation time (I/O included): 11 min, 49 sec, 191 millisec ************** 2.13 su macmini ? ale_loggamma(0,10007) Precomputation loggamma-values for q = 10007 and saving on file Precomputation time (I/O included): 0 min, 0 sec, 128 millisec ? preclngamma(0,10007) Precomputation lngamma-values for q = 10007 and saving on file Precomputation time (I/O included): 0 min, 0 sec, 212 millisec -------- ? ale_loggamma(0,305741) Precomputation loggamma-values for q = 305741 and saving on file Precomputation time (I/O included): 0 min, 3 sec, 847 millisec ? preclngamma(0,305741) Precomputation lngamma-values for q = 305741 and saving on file Precomputation time (I/O included): 0 min, 6 sec, 508 millisec -------- ? ale_loggamma(0,6766811) Precomputation loggamma-values for q = 6766811 and saving on file Precomputation time (I/O included): 1 min, 26 sec, 234 millisec ? preclngamma(0,6766811) Precomputation lngamma-values for q = 6766811 and saving on file Precomputation time (I/O included): 2 min, 24 sec, 972 millisec -------- ? ale_loggamma(0,10000019) Precomputation loggamma-values for q = 10000019 and saving on file Precomputation time (I/O included): 2 min, 8 sec, 653 millisec ? preclngamma(0,10000019) Precomputation lngamma-values for q = 10000019 and saving on file Precomputation time (I/O included): 3 min, 33 sec, 573 millisec -------- ? ale_loggamma(0,28227761) Precomputation loggamma-values for q = 28227761 and saving on file Precomputation time (I/O included): 6 min, 3 sec, 120 millisec ? preclngamma(0,28227761) Precomputation lngamma-values for q = 28227761 and saving on file Precomputation time (I/O included): 10 min, 1 sec, 711 millisec *************** 2.11.4 su macmini gp2c-run -pmy_ -g -W ale_loggamma.gp Questo script e' piu' veloce di log(Gamma) implementato in parigp ? ale_loggamma(0,10007) Precomputation loggamma-values for q = 10007 and saving on file Precomputation time (I/O included): 0 min, 0 sec, 136 millisec ? preclngamma(0,10007) Precomputation lngamma-values for q = 10007 and saving on file Precomputation time (I/O included): 0 min, 0 sec, 217 millisec --------- ? ale_loggamma(0,305741) Precomputation loggamma-values for q = 305741 and saving on file Precomputation time (I/O included): 0 min, 4 sec, 15 millisec ? preclngamma(0,305741) Precomputation lngamma-values for q = 305741 and saving on file Precomputation time (I/O included): 0 min, 6 sec, 542 millisec ---- ? ale_loggamma(0,6766811) Precomputation loggamma-values for q = 6766811 and saving on file Precomputation time (I/O included): 1 min, 31 sec, 460 millisec ? preclngamma(0,6766811) Precomputation lngamma-values for q = 6766811 and saving on file Precomputation time (I/O included): 2 min, 28 sec, 379 millisec ---- ? ale_loggamma(0,10000019) Precomputation loggamma-values for q = 10000019 and saving on file Precomputation time (I/O included): 2 min, 13 sec, 907 millisec ? preclngamma(0,10000019) Precomputation lngamma-values for q = 10000019 and saving on file Precomputation time (I/O included): 3 min, 36 sec, 43 millisec ------- ? ale_loggamma(0,28227761) Precomputation loggamma-values for q = 28227761 and saving on file Precomputation time (I/O included): 6 min, 14 sec, 590 millisec ? preclngamma(0,28227761) Precomputation lngamma-values for q = 28227761 and saving on file Precomputation time (I/O included): 10 min, 5 sec, 716 millisec *******/