mersenneforum.org

mersenneforum.org (https://www.mersenneforum.org/index.php)
-   Math (https://www.mersenneforum.org/forumdisplay.php?f=8)
-   -   pari's capabality (https://www.mersenneforum.org/showthread.php?t=16930)

devarajkandadai 2012-06-23 05:05

pari's capabality
 
I tried running the following program

{p(n)=isprime(n^2+1))}

It runs satisfactorily only till about n= 5000; beyond that pari does not seem to be able to handle.

science_man_88 2012-06-23 11:46

[QUOTE=devarajkandadai;303066]I tried running the following program

{p(n)=isprime(n^2+1))}

It runs satisfactorily only till about n= 5000; beyond that pari does not seem to be able to handle.[/QUOTE]

what version might be useful to some others here ( not me personally)? both versions I have seem to work for up to n=50000 at least.

science_man_88 2012-06-23 15:03

[QUOTE=science_man_88;303081]what version might be useful to some others here ( not me personally)? both versions I have seem to work for up to n=50000 at least.[/QUOTE]

doh I see now it works just slower than you wanted I guess. one thing I see that could speed speed it up is the fact that x^2 mod 6 -> 1,4,3,4,1,0 repeating so for x^2+1 to be 1 or 5 mod 6 ( the only ones that can be prime for numbers >3) if(x mod 2 == 1,return 0)

mart_r 2012-06-24 13:13

1 Attachment(s)
It works fine for me.

science_man_88 2012-06-24 17:34

[QUOTE=devarajkandadai;303066]I tried running the following program

{p(n)=isprime(n^2+1))}

It runs satisfactorily only till about n= 5000; beyond that pari does not seem to be able to handle.[/QUOTE]

I'm guessing you've already picked up on the fact that n<n^2+1 which sounds trivial but can be used, this shows that when y^2+1 is prime all indexes greater than that that fall in the groups y or -y mod y^2+1 are divisible by it I believe so as you find primes you can eliminate a lot from the search without primes outside the sequence.

science_man_88 2012-06-24 18:43

I've got a code that runs in 25-27 seconds in 2.4.2 and 19-21 seconds in 2.5.1 by the looks of it:

[CODE]a=vector(9000000,x,1);for(y=1,#a,if(a[y]==1,if(isprime(y^2+1),forstep(z=(y^2+1)-y,#a,[2*y,(y^2+1)-(2*y)],if(z%(y^2+1)==y || z%(y^2+1)==-y,a[z]=0)),a[y]=0),next()));a[/CODE]

only one thing to add to get it to print something other than 1's and 0's and it can add about 3 seconds to the run time. okay maybe I tested the add-on with 50000 still under 4 minutes in 2.5.1 doesn't sound bad for doing a search and print in a 9000000 element vector.


All times are UTC. The time now is 17:07.

Powered by vBulletin® Version 3.8.11
Copyright ©2000 - 2021, Jelsoft Enterprises Ltd.