Skip to content

Instantly share code, notes, and snippets.

@Hermann-SW
Created March 23, 2024 11:11
Show Gist options
  • Save Hermann-SW/b373c166d3a9e35253f3f4eab020b002 to your computer and use it in GitHub Desktop.
Save Hermann-SW/b373c166d3a9e35253f3f4eab020b002 to your computer and use it in GitHub Desktop.
Analyze Proth prine k*2^m+1 (odd k < 2^m, m≥1) distribution of "k"s for a given 1≤m≤20
P=extern("curl -s https://oeis.org/A080076/a080076.json.txt");
p2=eval(getenv("p2"));
r=eval(getenv("r"));
doit(r)={
for(m=r[1],r[2],
cnt=Vec(0,2^p2);
foreach(P,p,if(2^m<p&&p<2^(2*m),cnt[1+(p>>m)>>(m-p2)]++));
print(m," ",cnt," (0..",2^(p2)-1,") vecsum=",vecsum(cnt));
)
};
doit(r);
@Hermann-SW
Copy link
Author

Frequency of Proth prime "k"s for a given m much bigger on the leftmost part (here for 2^p2==8 parts):

pi@raspberrypi5:~ $ p2=3 r=[15,20] gp -q < proth_m.gp 
15 [1937, 850, 397, 402, 411, 391, 394, 387] (0..7) vecsum=5169
16 [3564, 1583, 780, 750, 733, 784, 756, 747] (0..7) vecsum=9697
17 [6651, 3020, 1436, 1414, 1440, 1474, 1403, 1345] (0..7) vecsum=18183
18 [12483, 5662, 2733, 2744, 2720, 2617, 2682, 2618] (0..7) vecsum=34259
19 [23588, 10637, 5195, 5102, 5055, 5108, 5030, 5014] (0..7) vecsum=64729
20 [44451, 20207, 9921, 9902, 9605, 9596, 9412, 9395] (0..7) vecsum=122489
pi@raspberrypi5:~ $ 

@Hermann-SW
Copy link
Author

With C set to this array of length 2^p2==1024 for m=20 ...

pi@raspberrypi5:~ $ p2=10 r=[20,20] gp -q < proth_m.gp 
20 [4978, 1530, 1517, ... , 88, 75, 74, 60] (0..1023) vecsum=122489
pi@raspberrypi5:~ $ 

... this command draws diagram to reveal more of the Proth prime k frequency distribution for m=20:

? plotinit(0,#C,800);plotscale(0,0,#C,0,C[1]);plotlines(0,[1..#C],[C[X]|X<-[1..#C]]);plotlines(0,[0,#C],[0,0]);for(h=1,4,plotlines(0,[2^(10-2*h),2^(12-2*h)],[10*h,10*h]));plotdraw(0)
? 

Proth pari_plot
The bottom line is x-axis for 1024 buckets, y-axis is number of "k"s falling into each bucket.
The other horizontal straight lines indicate 4^i steps of frequency.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment