Re: Big GOTCHAs WITH forprime()

[Bill Allombert <Bill.Allombert@math.u-bordeaux.fr>]

On Thu, Jan 25, 2024 at 04:40:07AM -0800, Ilya Zakharevich wrote:
> On Thu, Jan 25, 2024 at 10:39:35AM +0100, Bill Allombert wrote:
> > >   https://pari.math.u-bordeaux.fr/cgi-bin/bugreport.cgi?bug=2520#10
> > 
> > Yes, I was running the immediate running version but without
> > primelimit=400000000000
> > 
> > Now on linux, with GP 2.15.4:
> > 
> > % gp -qf --default primelimit=400000000000
> > ? my((f(s)=forprime(p=s,s+10^8,)),n=98,x=190333*10^12,y=x+2*10^12,t,z,t0);while(n>0,t=getwalltime();n-=1+if(if(1,f(z=(x+3*y)/4);20000<(t0=-(t-(t=getwalltime())))),x=z;4,y=z;0);print([1.*x,1.*y,t0,n]))
> > 
> > [190333000000000000.00000000000000000000, 190334500000000000.00000000000000000000, 10650, 97]
> > [190334125000000000.00000000000000000000, 190334500000000000.00000000000000000000, 32178, 92]
> > [190334125000000000.00000000000000000000, 190334406250000000.00000000000000000000, 10259, 91]
> > [190334125000000000.00000000000000000000, 190334335937500000.00000000000000000000, 10436, 90]
> > [190334125000000000.00000000000000000000, 190334283203125000.00000000000000000000, 10284, 89]
> > [190334125000000000.00000000000000000000, 190334243652343750.00000000000000000000, 10629, 88]
> > [190334125000000000.00000000000000000000, 190334213989257812.50000000000000000000, 10562, 87]
> 
> Thanks for reporting back!  However:
> 
> I repeat the 4th time:
>      This is not the code mentioned above, in 2520#10.

Ah sorry! For some reason, I kept pasting the wrong line.

So it seems the first bad commit is

commit a286058de199e0477bb51d82a3f092ca71a47ad9
Author: Karim Belabas <Karim.Belabas@math.u-bordeaux1.fr>
Date:   Mon May 26 11:56:00 2014 +0200

    25- (gp -p N) or (primelimit=N in gprc_ for N >= 436273290 resulted in an
       incorrect primetable. N.B. Such commands are now useless: needed primes
       are produced dynamically anyway.


gdb --arg ./gp.dbg -p400G
parisizemax = 4000002048, primelimit = 400000000000
? my((f(s)=forprime(p=s,s+10^8,)),n=67,x=190334138131456478.73170673847198486328,y=190334138444662596,t,z,t0);while(n>0,t=getabstime();n-=1+if(if(1,f(z=(x+3*y)/4);20000<(t0=-(t-(t=getabstime())))),x=z;4,y=z;0);print([1.*x,1.*y,t0,n]))

program received signal SIGSEGV, Segmentation fault.
0x0000555555baac69 in sieve_block (a=190334138381396059, b=190334138381914081, maxpos=32376,
    sieve=0x7ffff6f200f8 "\377\367\337\377\377\377\377\377\337\377\377\373\377\377\377\377\277\377\367~\375\377\377\333\377\376\377\377\277\377\377\377\375\377\377\377\377\277\377\377\377\377\377\377\357\373\377\377\367\377\337\377\377\377\377\376\373\377~\373\377\377\377\376\377\377\367\357\377\373\277\377\377\377\377\277\377\375\377\377\377\377\377\375\376\375\377\377\376\377\372\357\377\273\377\337\367\377\377\367\277\377\377\377\377\377\377\337\377\377\377\377\357\377\377\277\371\377\377\377\376\377\335\377g\377\377\377\377\367\177\377\377\377\377\377\377\377\377\377}\377\357\377\377\377\377\337\377\377\377\177\375\377\377\377\377\377\377\377\377\377\367\366\357\377\177\375\377\367\377\377\377\373\377=\377\377\377\377\377\377\377\377\377\315\377\377\377\377\377\377\377\377\377\277\377\377\377", <incomplete
sequence \357>...) at ../src/language/forprime.c:781
781         NEXT_PRIME_VIADIFF(p, d); /* starts at p = 3 */
(gdb) p p
$1 = 436273009
(gdb) p d
$2 = (byteptr) 0x7fff087c2001 <error: Cannot access memory at address 0x7fff087c2001>

the "obvious" fix is

-    NEXT_PRIME_VIADIFF(p, d); /* starts at p = 3 */
+    NEXT_PRIME_VIADIFF_CHECK(p, d); /* starts at p = 3 */

which leads to

? my((f(s)=forprime(p=s,s+10^8,)),n=67,x=190334138131456478.73170673847198486328,y=190334138444662596,t,z,t0);while(n>0,t=getabstime();n-=1+if(if(1,f(z=(x+3*y)/4);20000<(t0=-(t-(t=getabstime())))),x=z;4,y=z;0);print([1.*x,1.*y,t0,n]))

  ***   at top-level: ...etabstime();n-=1+if(if(1,f(z=(x+3*y)/4);20000<
  ***                                             ^---------------------
  ***   in function f: forprime(p=s,s+10^8,)
  ***                  ^---------------------
  ***   not enough precomputed primes

which is not quite accurate.

Cheers,
Bill.