How do you efficiently sieve for prime 3/4-tuples?

[Puzzle-Peter]

Originally Posted by R. Gerbicz

I'm in the process to rewrite gsieve.c for your new problem. It will be possible to sieve large primes also (so q>2^31), what gsieve currently doesn't handle.

Wow, this is really great news! Thanks so much for the effort you're putting into this! I will of course take a close look at the new code and hopefully learn a few things from doing so.

Again, thank you!

[henryzz]

Originally Posted by R. Gerbicz

I'm in the process to rewrite gsieve.c for your new problem. It will be possible to sieve large primes also (so q>2^31), what gsieve currently doesn't handle.

This does sound good. Can I suggest that comments in the code might help others understand and help improve the code in the future?
One improvement I noticed(untested) is that on the time critical lines:
B[h>>6]&=~Bits[h&63];
A[b>>6]&=~Bits[b&63];
You calculate the complement ~ each time. Since this is such time critical code(from memory most of the runtime between these two lines) would it make sense to store in Bits the complement of the current values and on other lines using Bits(not many and I think not time critical) calculate the complement. If necessary you could have two Bits arrays, one containing the complement.
I don't know how much this would speed things up but it might take some significant time off some jobs.

[R. Gerbicz]

I have just completed my new program, see: https://sites.google.com/site/robertgerbicz/polysieve
to use it for your problem, use this:

Give the degree of the P polynom: 1
Give the 0-th coefficient of P: -4
Give the 1-th coefficient of P: -2
Give the degree of the Q polynom: 3
Give the 0-th coefficient of Q: 0
Give the 1-th coefficient of Q: 0
Give the 2-th coefficient of Q: 2
Give the 3-th coefficient of Q: 1
Give the number of c values for the sieve: 3
0-th c value: 1
1-th c value: -1
2-th c value: 5

As promised it can sieve large primes, really fast. For larger range it is faster up to 2310*10^9, in the best setup use such range that is a multiple of 2310*10^9, or you will loose roughly one hour of work or so depending on your computer. The memory usage is somewhat near to gsieve on small primes. For larger primes it will use even less memory, but this depends on the size of hashtable (the default size is 256MB, you can increase this).

You can see that it is a more general sieving code, and can handle any polynomial problem of this form. Ask if something is not clear.

Just to make clear: on your problem you need to find s, where s,s+2 is a twin prime. But the sieve needs s=k*b^n+d, so to find an easy primality proof also on s,s+2 it means that d=-1. One such example is s=2673*2^116-1. For the code it is a little faster if you use b=2. As this is only a sieving code, it does NOT check if s and s+2 is prime or not, that is your task.
Use k, where k<2^63.

ps. Precomputed the ~Bits[] array, thanks Henryzz.

[Puzzle-Peter]

That was fast work! I'll give it a try and post feedback & timings as soon as possible.

As I said I plan to use already known pairs of twins s, s+2 and the vast majority of these is of the form k*2^n +-1, so that fits perfectly.

Thanks a lot!

[R. Gerbicz]

Originally Posted by Puzzle-Peter

That was fast work! I'll give it a try and post feedback & timings as soon as possible.

Too fast, discovered a bug, for some prime divisors of Q(s) the sieve is broken. However this has no effect on our current polynom. I will correct this soon,

[henryzz]

Been fiddling today. I am very impressed. Ran a long test where the second sieve took much longer than expected to begin with but sped up as it sieved deeper. It especially seemed to speed up after the hash updates. Is there a particular reason for how the hash updates are distributed or could it be done more often? Obviously too often would be harmful as it does take a little time. I am going to rerun this test tomorrow with a larger hashtable to see if that improves things.

Here is the log of the slow run:

Sieve P(s)+a*Q(s)+c for multiple c values, with fixed s=k*b^n+d; P,Q is polynom.
Give k: 467343
Give b: 2
Give n: 1000
Give d: -1
Give the degree of the P polynom: 1
Give the 0-th coefficient of P: -4
Give the 1-th coefficient of P: -2
Give the degree of the Q polynom: 3
Give the 0-th coefficient of Q: 0
Give the 1-th coefficient of Q: 0
Give the 2-th coefficient of Q: 2
Give the 3-th coefficient of Q: 1
Give the number of c values for the sieve: 3
0-th c value: 1
1-th c value: -1
2-th c value: 5
Give start and end value for 'a' (in billions)! 0
2310
Give the limit for sieving primes (maxp): 10000000000
Give the name of the file to output the numbers! 467343.pfgw
Using primes for wheelsieve up to 11
initialization of all small primes done, primepi(1000000000)=50847534
initialization of residue and inverse array done up to p=1000000000.
In wheel sieve there will be 64 different combinations mod 2310
Start sieving for a=0..2310000000000.
64/64; #candidates=131315557; 4628 sec. so far.
Sort the 'a' array.
Done the sorting. Sieve for large primes.
Use M=6;mres=2;
Sieved up to 2000000000; #candidates=118975467; time so far=13468 sec.
Update the hash table...done
Sieved up to 3000000000; #candidates=112465624; time so far=18088 sec.
Sieved up to 4000000000; #candidates=108132559; time so far=21820 sec.
Update the hash table...done
Sieved up to 5000000000; #candidates=104925658; time so far=24372 sec.
Sieved up to 6000000000; #candidates=102395988; time so far=26425 sec.
Sieved up to 7000000000; #candidates=100322524; time so far=28192 sec.
Sieved up to 8000000000; #candidates=98572783; time so far=29901 sec.
Update the hash table...done
Sieved up to 9000000000; #candidates=97060698; time so far=31049 sec.
Sieved up to 10000000000; #candidates=95733174; time so far=32077 sec.
Output the remaining 'a' values to the abc file, wait...
Done. Number of remaining candidates=95733174.
Completed in 32146 sec.

Here is a run 1/11th of the size. It runs in 1/16.45th the time even though the first sieve runs at 8/11 the speed.

Sieve P(s)+a*Q(s)+c for multiple c values, with fixed s=k*b^n+d; P,Q is polynom.
Give k: 467343
Give b: 2
Give n: 1000
Give d: -1
Give the degree of the P polynom: 1
Give the 0-th coefficient of P: -4
Give the 1-th coefficient of P: -2
Give the degree of the Q polynom: 3
Give the 0-th coefficient of Q: 0
Give the 1-th coefficient of Q: 0
Give the 2-th coefficient of Q: 2
Give the 3-th coefficient of Q: 1
Give the number of c values for the sieve: 3
0-th c value: 1
1-th c value: -1
2-th c value: 5
Give start and end value for 'a' (in billions)! 0
210
Give the limit for sieving primes (maxp): 10000000000
Give the name of the file to output the numbers! 467343_7_d12.pfgw
Using primes for wheelsieve up to 7
initialization of all small primes done, primepi(1000000000)=50847534
initialization of residue and inverse array done up to p=1000000000.
In wheel sieve there will be 8 different combinations mod 210
Start sieving for a=0..210000000000.
8/8; #candidates=11934182; 520 sec. so far.
Sort the 'a' array.
Done the sorting. Sieve for large primes.
Use M=6;mres=2;
Sieved up to 2000000000; #candidates=10814744; time so far=737 sec.
Update the hash table...done
Sieved up to 3000000000; #candidates=10222163; time so far=891 sec.
Sieved up to 4000000000; #candidates=9828206; time so far=1021 sec.
Update the hash table...done
Sieved up to 5000000000; #candidates=9535643; time so far=1140 sec.
Sieved up to 6000000000; #candidates=9305757; time so far=1249 sec.
Sieved up to 7000000000; #candidates=9117937; time so far=1354 sec.
Sieved up to 8000000000; #candidates=8958570; time so far=1454 sec.
Update the hash table...done
Sieved up to 9000000000; #candidates=8820961; time so far=1551 sec.
Sieved up to 10000000000; #candidates=8700451; time so far=1646 sec.
Output the remaining 'a' values to the abc file, wait...
Done. Number of remaining candidates=8700451.
Completed in 1652 sec.

I also found my first triple using this:
-2*(467343*2^1000-1)-4+1788349880*((467343*2^1000-1)^3+2*(467343*2^1000-1)^2)+1
-2*(467343*2^1000-1)-4+1788349880*((467343*2^1000-1)^3+2*(467343*2^1000-1)^2)-1
-2*(467343*2^1000-1)-4+1788349880*((467343*2^1000-1)^3+2*(467343*2^1000-1)^2)+5
I will prove it tomorrow.

I suppos with P(s)=0 and Q(s)=2^n this would effectively sieve the same form as gsieve. Would/should it complain if s isn't used? The other question is how much overhead is there from supporting the wider forms. I will hopefully test this tommorrow and compare the results of gsieve with this(Is this a bugged Q(s) as mentioned above?)

[R. Gerbicz]

Originally Posted by henryzz

Been fiddling today. I am very impressed. Ran a long test where the second sieve took much longer than expected to begin with but sped up as it sieved deeper. It especially seemed to speed up after the hash updates. Is there a particular reason for how the hash updates are distributed or could it be done more often? Obviously too often would be harmful as it does take a little time. I am going to rerun this test tomorrow with a larger hashtable to see if that improves things.

The updates comes in geometric progression with quotient=2. I will add an option to modify the frequency. Larger hashtable helps, though even with 256MB of table and 100 million survivors that means you will need a binary search with roughly only 1/21 probability,

Originally Posted by henryzz

I suppos with P(s)=0 and Q(s)=2^n this would effectively sieve the same form as gsieve. Would/should it complain if s isn't used? The other question is how much overhead is there from supporting the wider forms. I will hopefully test this tommorrow and compare the results of gsieve with this(Is this a bugged Q(s) as mentioned above?)

Note that here the coefficients should be 64 signed bits, so choosing 2^n for n>=64 is impossible. In wider forms there is no much overhead, the overhead depends only on the degree (what should be less than 100).

[henryzz]

Originally Posted by R. Gerbicz

Note that here the coefficients should be 64 signed bits, so choosing 2^n for n>=64 is impossible. In wider forms there is no much overhead, the overhead depends only on the degree (what should be less than 100).

I knew I had a use for s lol. s=1*2^n-0 should work I think?

[R. Gerbicz]

Originally Posted by henryzz

I knew I had a use for s lol. s=1*2^n-0 should work I think?

Oh, yes. s=1*2^n-0 will work, with P(s)=0, so degree=0 and coefficient=0 for P.

[Puzzle-Peter]

Couldn't do a lot of tests yet, but I have one result. Using "my" polynome and an s of ~11K decimal digits, I first sieved a=0 to 1G up to pmax=10G using the original parameters. Running time was 4713s. I then reduced SmallPrimes from 1G to 100M and got a running time of 864s, the output files were identical.

I'll try a=0 to 2310G next.

[henryzz]

If you want smaller efficient ranges to sieve then reduce bound_small_primes to 7 and you can sieve in multiples of 210G. bound_small_primes at 5 sieves in multiples of 30G.

For a range of 2310G I think that it is more efficient to increase the hashtable size from 31 bits. 32 is running much faster today than 31 yesterday. I will need to check that the things I had running yesterday weren't affecting it much. It looks like it will finish in <23k seconds