mersenneforum.org  

Go Back   mersenneforum.org > Great Internet Mersenne Prime Search > Math

Reply
 
Thread Tools
Old 2018-03-17, 11:44   #45
science_man_88
 
science_man_88's Avatar
 
"Forget I exist"
Jul 2009
Dumbassville

8,369 Posts
Default

Quote:
Originally Posted by Dr Sardonicus View Post
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.

Last fiddled with by science_man_88 on 2018-03-17 at 12:25
science_man_88 is offline   Reply With Quote
Old 2018-03-17, 15:14   #46
Dr Sardonicus
 
Dr Sardonicus's Avatar
 
Feb 2017
Nowhere

3×7×132 Posts
Default

Quote:
Originally Posted by Dr Sardonicus View Post
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.

Last fiddled with by Dr Sardonicus on 2018-03-17 at 15:23
Dr Sardonicus is offline   Reply With Quote
Old 2018-03-17, 17:00   #47
science_man_88
 
science_man_88's Avatar
 
"Forget I exist"
Jul 2009
Dumbassville

8,369 Posts
Default

Quote:
Originally Posted by Dr Sardonicus View Post
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.

Last fiddled with by science_man_88 on 2018-03-17 at 17:02
science_man_88 is offline   Reply With Quote
Old 2018-03-18, 15:55   #48
Dr Sardonicus
 
Dr Sardonicus's Avatar
 
Feb 2017
Nowhere

3×7×132 Posts
Default

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
Dr Sardonicus is offline   Reply With Quote
Old 2018-03-18, 16:26   #49
science_man_88
 
science_man_88's Avatar
 
"Forget I exist"
Jul 2009
Dumbassville

8,369 Posts
Default

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.
science_man_88 is offline   Reply With Quote
Old 2018-03-18, 17:12   #50
Batalov
 
Batalov's Avatar
 
"Serge"
Mar 2008
Phi(4,2^7658614+1)/2

3·3,041 Posts
Default

Just in case, if anyone is sieiving or testing for m=17, a(17) > 238.
Batalov is offline   Reply With Quote
Old 2018-03-18, 17:20   #51
paulunderwood
 
paulunderwood's Avatar
 
Sep 2002
Database er0rr

22·859 Posts
Default

Quote:
Originally Posted by Batalov View Post
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
paulunderwood is offline   Reply With Quote
Old 2018-03-18, 19:47   #52
Dr Sardonicus
 
Dr Sardonicus's Avatar
 
Feb 2017
Nowhere

3·7·132 Posts
Default

Quote:
Originally Posted by science_man_88 View Post
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...

Last fiddled with by Dr Sardonicus on 2018-03-18 at 19:53
Dr Sardonicus is offline   Reply With Quote
Old 2018-03-18, 20:04   #53
a1call
 
a1call's Avatar
 
"Rashid Naimi"
Oct 2015
Remote to Here/There

19×101 Posts
Default

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
a1call is offline   Reply With Quote
Old 2018-03-18, 20:19   #54
science_man_88
 
science_man_88's Avatar
 
"Forget I exist"
Jul 2009
Dumbassville

8,369 Posts
Default

Quote:
Originally Posted by Dr Sardonicus View Post
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.

Last fiddled with by science_man_88 on 2018-03-18 at 20:36
science_man_88 is offline   Reply With Quote
Old 2018-03-18, 20:47   #55
a1call
 
a1call's Avatar
 
"Rashid Naimi"
Oct 2015
Remote to Here/There

19·101 Posts
Default

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.
a1call is offline   Reply With Quote
Reply

Thread Tools


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 k*a(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

All times are UTC. The time now is 01:08.

Tue Oct 20 01:08:36 UTC 2020 up 39 days, 22:19, 0 users, load averages: 1.89, 2.03, 2.13

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

This forum has received and complied with 0 (zero) government requests for information.

Permission is granted to copy, distribute and/or modify this document under the terms of the GNU Free Documentation License, Version 1.2 or any later version published by the Free Software Foundation.
A copy of the license is included in the FAQ.