Smallest prime of the form a^2^m + b^2^m, m>=14 - Page 5

science_man_88 · 2018-03-17, 11:44

Quote:

Originally Posted by Dr Sardonicus

The number of pairs (a, b) with 1 <= a <= b <= N, and gcd(a, b) = 1, is SUM, b = 1 to N, eulerphi(b). This is the number of fractions in the Farey sequence of order N.

Asymptotically, this is 6/pi^2 * N^2. Long known.

The additional condition that a + b is odd, knocks out the pairs (a, b) with a and b both odd. In this case, b - a is even and gcd(b - a, b) = 1, so for odd b > 1, the condition knocks out exactly 1/2*eulerphi(b) pairs. For b = 1 it knocks out the pair (1, 1).

I'm too lazy to work out the modified asymptotic. I'm sure it is long known, but alas I'm also too lazy to look it up
;-)

Using me and serges calculations it seems to be roughly 0.2*N^2 I get a>11 or b>10 or both for the next one.

Dr Sardonicus · 2018-03-17, 15:14

Quote:

Originally Posted by Dr Sardonicus

The number of pairs (a, b) with 1 <= a <= b <= N, and gcd(a, b) = 1, is SUM, b = 1 to N, eulerphi(b). This is the number of fractions in the Farey sequence of order N.

Asymptotically, this is 6/pi^2 * N^2. Long known.

The additional condition that a + b is odd, knocks out the pairs (a, b) with a and b both odd. In this case, b - a is even and gcd(b - a, b) = 1, so for odd b > 1, the condition knocks out exactly 1/2*eulerphi(b) pairs. For b = 1 it knocks out the pair (1, 1).

I'm too lazy to work out the modified asymptotic. I'm sure it is long known, but alas I'm also too lazy to look it up
;-)

I screwed up. Mea culpa.

The coefficient 6/pi^2 is for the number of pairs of positive integer (a, b), both less than or equal to N, with gcd(a, b) = 1 -- but without the requirement a < b. This coefficient may be viewed as the limiting value, as N increases without bound, of the probability that two positive integers, chosen at random in [1, N], are relatively prime. By considering separately the possible prime common factors 2, 3, 5, etc and a bit of handwaving, you obtain the coefficient as being the product

(1 - 1/2^2)*(1 - 1/3^2)*(1 - 1/5^2)...

whose inverse is obviously the sum, from k = 1 to infinity, of 1/k^2, whose value is well known to be pi^2/6. Thus, the coefficient is 6/pi^2.

Imposing the requirement a < b obviously cuts the number of pairs in half, giving the correct asymptotic for the Farey sequence of order N as 3/pi^2 * N^2. The coefficient 3/pi^2 is

0.30396355092701331433163838962918291671, approximately.

Being too lazy either to work out or look up the variant with the further requirement that a and b be of different parity, I simply had Pari-GP find the actual value numerically up to some ridiculous limit, and print out the value every so often. It was interesting to watch as an increasing number of decimal digits remained fixed, but the convergence was slow, so I Ctrl-C'd it out of my misery.

The asymptotic for the variant clearly is 2/pi^2 * N^2, the coefficient 2/pi^2 being

0.20264236728467554288775892641945527781 approximately.

science_man_88 · 2018-03-17, 17:00

Quote:

Originally Posted by Dr Sardonicus

I screwed up. Mea culpa.

The coefficient 6/pi^2 is for the number of pairs of positive integer (a, b), both less than or equal to N, with gcd(a, b) = 1 -- but without the requirement a < b. This coefficient may be viewed as the limiting value, as N increases without bound, of the probability that two positive integers, chosen at random in [1, N], are relatively prime. By considering separately the possible prime common factors 2, 3, 5, etc and a bit of handwaving, you obtain the coefficient as being the product

(1 - 1/2^2)*(1 - 1/3^2)*(1 - 1/5^2)...

whose inverse is obviously the sum, from k = 1 to infinity, of 1/k^2, whose value is well known to be pi^2/6. Thus, the coefficient is 6/pi^2.

Imposing the requirement a < b obviously cuts the number of pairs in half, giving the correct asymptotic for the Farey sequence of order N as 3/pi^2 * N^2. The coefficient 3/pi^2 is

0.30396355092701331433163838962918291671, approximately.

Being too lazy either to work out or look up the variant with the further requirement that a and b be of different parity, I simply had Pari-GP find the actual value numerically up to some ridiculous limit, and print out the value every so often. It was interesting to watch as an increasing number of decimal digits remained fixed, but the convergence was slow, so I Ctrl-C'd it out of my misery.

The asymptotic for the variant clearly is 2/pi^2 * N^2, the coefficient 2/pi^2 being

0.20264236728467554288775892641945527781 approximately.

To be fair, my "calculation" was done with PARI/gp code that counted them for me. I still wish I could figure my additive inverse version out. If it weren't for the 14,15 cases it seemed most often a composite m had a odd and b even. Of course so did 7( and 8 breaks the general composite rule). Just looking at the data.

Dr Sardonicus · 2018-03-18, 15:55

Just for practice, I (tried to) write Pari-GP scripts to perform various tasks associated with this problem. I compiled a vector of all pairs [a, b] with b up to 200 for which

a < b, gcd(a, b) = 1, and a + b is odd.

I deliberately left the testing of the conditions appallingly clumsy. Since I didn't bother to compute in advance how many pairs there were going to be, I used concat() to append each new pair to my vector. That took way more time than computing the pairs [a, b]. No matter, I had my vector of pairs, 8156 of them. If I had done the computation of the number of pairs in advance, this part of the operation would have been much faster.

Next, I decided to see how effective sieving might be. Using the original exponent e = 2^14, I compiled a vector of primes p < 2^26, p congruent to 1 (mod 2^15). There are 228 of them.

I then wrote a script that initialized a count at 0 and went through my vector of pairs [a,b]. If my script is writ right (I invite anyone to check the results, I'm not good at writing scripts), for each pair [a,b] it goes through my vector of primes p, and computes

r = Mod(a,p)^e + Mod(b,p)^e.

Then if r is zero, it increments the count, and immediately moves on to the next pair [a,b].

(I note that the smallest value 2^(2^14) + 1 of a^e + b^e is way bigger than my prime limit 2^26.)

At the end of the run, the count of pairs sieved out was 4065. Now 4065./8156 = 0.4984, approximately. The "grasping at straws" product of (1. - 2^14/p), taken over the 228 primes on my list, is 0.4995, approximately. That is much closer agreement than I had hoped for!

science_man_88 · 2018-03-18, 16:26

my code for counting them was just:

Code:

my(a=0);for(x=1,200,forstep(y=x-1,1,-2,if(gcd(x,y)==1,a++)));a

Done in about 50 ms I guess. Where as the version using concat was 19591 ms I may be able to do it with two vectors with fewer concats to the main vector though.

Batalov · 2018-03-18, 17:12

Just in case, if anyone is sieiving or testing for m=17, a(17) > 238.

paulunderwood · 2018-03-18, 17:20

Quote:

Originally Posted by Batalov

Just in case, if anyone is sieiving or testing for m=17, a(17) > 238.

That's 311503 decimal digits -- PRP'ing with generic mod reduction is really going to bite

Dr Sardonicus · 2018-03-18, 19:47

Quote:

Originally Posted by science_man_88

my code for counting them was just:

Code:

my(a=0);for(x=1,200,forstep(y=x-1,1,-2,if(gcd(x,y)==1,a++)));a

Done in about 50 ms I guess. Where as the version using concat was 19591 ms I may be able to do it with two vectors with fewer concats to the main vector though.

I went back and did the counting in advance (up to the limit b = 200) thusly:

Code:

? s=0;r=1;for(i=2,200,if(r,t=eulerphi(i),t=eulerphi(i)/2);s+=t;r=!r);print(s)
8156
time = 0 ms.

For creating the vector of pairs [a, b] I then did this: I'm sure the testing of whether i+j is odd can be handled more expeditiously, but it probably wouldn't make much difference in the time.

Code:

? k=1;v=vector(s);for(j=2,200,for(i=1,j-1,if((i+j)%2==1&&gcd(i,j)==1,v[k]=[i,j];k++)))
time = 47 ms.
? print(v[s])
[199, 200]

It looks to me like improving times in looking for pseudoprimes would depend on doing as much sieving as possible. This would mean extending the "factor base" of primes congruent to 1 (mod 2*e) way, way beyond my tiny limit of 2^26. Doubling the log of the prime limit would perhaps sieve out half the pairs left by sieving to the old limit...

a1call · 2018-03-18, 20:04

My 2¢ worth,

The fastest filtering is size based since the objective is to find the smallest PRP.
You can set a narrow test range based on size, assign each band to a different thread and advance upwards after each completion. The size filtering before any trial by division or PRP test seems most efficient to me.
It should avoid a lot of unnecessary lager PRP tests.

science_man_88 · 2018-03-18, 20:19

Quote:

Originally Posted by Dr Sardonicus

I went back and did the counting in advance (up to the limit b = 200) thusly:

Code:

? s=0;r=1;for(i=2,200,if(r,t=eulerphi(i),t=eulerphi(i)/2);s+=t;r=!r);print(s)
8156
time = 0 ms.

For creating the vector of pairs [a, b] I then did this: I'm sure the testing of whether i+j is odd can be handled more expeditiously, but it probably wouldn't make much difference in the time.

Code:

? k=1;v=vector(s);for(j=2,200,for(i=1,j-1,if((i+j)%2==1&&gcd(i,j)==1,v[k]=[i,j];k++)))
time = 47 ms.
? print(v[s])
[199, 200]

It looks to me like improving times in looking for pseudoprimes would depend on doing as much sieving as possible. This would mean extending the "factor base" of primes congruent to 1 (mod 2*e) way, way beyond my tiny limit of 2^26. Doubling the log of the prime limit would perhaps sieve out half the pairs left by sieving to the old limit...

You can also treat composite and prime a values seperately. Odd primes always have all opposite parity less than them, composites don't always except maybe powers of two. But again doubtful this decreases the time. EDIT: you actually have that all odd numbers have the powers of two less than them. So really the ones in need of findng, are only when 4n+2 and odd numbers crash.

a1call · 2018-03-18, 20:47

Chances are you can eliminate some 10's of very expensive PRP tests by performing a GCD test between each candidate :

GCD (a^q+b^q,(a+k)^+(b+k')^q)

For arbitrary values of k and k'.

Such a test can speed things up, but gets less and less effective as the exponent/PRP gets larger.
There is not much saving if any, in doing more than one such test.

2018-03-18, 15:55	#48
Dr Sardonicus Feb 2017 Nowhere 11043₈ Posts	Just for practice, I (tried to) write Pari-GP scripts to perform various tasks associated with this problem. I compiled a vector of all pairs [a, b] with b up to 200 for which a < b, gcd(a, b) = 1, and a + b is odd. I deliberately left the testing of the conditions appallingly clumsy. Since I didn't bother to compute in advance how many pairs there were going to be, I used concat() to append each new pair to my vector. That took way more time than computing the pairs [a, b]. No matter, I had my vector of pairs, 8156 of them. If I had done the computation of the number of pairs in advance, this part of the operation would have been much faster. Next, I decided to see how effective sieving might be. Using the original exponent e = 2^14, I compiled a vector of primes p < 2^26, p congruent to 1 (mod 2^15). There are 228 of them. I then wrote a script that initialized a count at 0 and went through my vector of pairs [a,b]. If my script is writ right (I invite anyone to check the results, I'm not good at writing scripts), for each pair [a,b] it goes through my vector of primes p, and computes r = Mod(a,p)^e + Mod(b,p)^e. Then if r is zero, it increments the count, and immediately moves on to the next pair [a,b]. (I note that the smallest value 2^(2^14) + 1 of a^e + b^e is way bigger than my prime limit 2^26.) At the end of the run, the count of pairs sieved out was 4065. Now 4065./8156 = 0.4984, approximately. The "grasping at straws" product of (1. - 2^14/p), taken over the 228 primes on my list, is 0.4995, approximately. That is much closer agreement than I had hoped for! Last fiddled with by Dr Sardonicus on 2018-03-18 at 16:00 Reason: correcting wrong words

2018-03-18, 16:26	#49
science_man_88 "Forget I exist" Jul 2009 Dumbassville 20C0₁₆ Posts	my code for counting them was just: Code: my(a=0);for(x=1,200,forstep(y=x-1,1,-2,if(gcd(x,y)==1,a++)));a Done in about 50 ms I guess. Where as the version using concat was 19591 ms I may be able to do it with two vectors with fewer concats to the main vector though.

2018-03-18, 20:04	#53
a1call "Rashid Naimi" Oct 2015 Remote to Here/There 807₁₆ Posts	My 2¢ worth, The fastest filtering is size based since the objective is to find the smallest PRP. You can set a narrow test range based on size, assign each band to a different thread and advance upwards after each completion. The size filtering before any trial by division or PRP test seems most efficient to me. It should avoid a lot of unnecessary lager PRP tests. Last fiddled with by a1call on 2018-03-18 at 20:08

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
Is there a prime of the form......	PawnProver44	Miscellaneous Math	9	2016-03-19 22:11
OEIS A071580: Smallest prime of the form ka(n-1)a(n-2)...a(1)+1	arbooker	And now for something completely different	14	2015-05-22 23:18
Smallest prime with a digit sum of 911	Stargate38	Puzzles	6	2014-09-29 14:18
Smallest floor of k for cullen prime	Citrix	Prime Cullen Prime	12	2007-04-26 19:52
Smallest ten-million-digit prime	Heck	Factoring	9	2004-10-28 11:34

2018-03-18, 17:12	#50
Batalov "Serge" Mar 2008 Phi(4,2^7658614+1)/2 3⁶·13 Posts	Just in case, if anyone is sieiving or testing for m=17, a(17) > 238.

2018-03-18, 20:47	#55
a1call "Rashid Naimi" Oct 2015 Remote to Here/There 3×5×137 Posts	Chances are you can eliminate some 10's of very expensive PRP tests by performing a GCD test between each candidate : GCD (a^q+b^q,(a+k)^+(b+k')^q) For arbitrary values of k and k'. Such a test can speed things up, but gets less and less effective as the exponent/PRP gets larger. There is not much saving if any, in doing more than one such test.