Frobenius numbers
- Task
a(n) = prime(n)*prime(n+1) - prime(n) - prime(n+1), where prime(n) < 10,000
C#
<lang csharp>using System.Collections.Generic; using System.Linq; using static System.Console; using static System.Math;
class Program { static void Main() {
int c = 0, d = 0, f, lim = 1000000, l2 = lim / 100; var Frob = PG.Primes((int)Sqrt(lim) + 1).ToArray();
for (int n = 0, m = 1; m < Frob.Length; n = m++) {
if ((f = Frob[n] * Frob[m] - Frob[n] - Frob[m]) < l2) d++;
Write("{0,7:n0} {1}", f , ++c % 10 == 0 ? "\n" : ""); }
WriteLine("\n\nFound {0} Frobenius numbers of consecutive primes under {1:n0}," +
"of which {2} were under {3:n0}", c, lim, d, l2); } }
class PG { public static IEnumerable<int> Primes(int lim) {
var flags = new bool[lim + 1]; int j; yield return 2;
for (j = 4; j <= lim; j += 2) flags[j] = true; j = 3;
for (int d = 8, sq = 9; sq <= lim; j += 2, sq += d += 8)
if (!flags[j]) { yield return j;
for (int k = sq, i = j << 1; k <= lim; k += i) flags[k] = true; }
for (; j <= lim; j += 2) if (!flags[j]) yield return j; } }</lang>
- Output:
1 7 23 59 119 191 287 395 615 839 1,079 1,439 1,679 1,931 2,391 3,015 3,479 3,959 4,619 5,039 5,615 6,395 7,215 8,447 9,599 10,199 10,811 11,447 12,095 14,111 16,379 17,679 18,767 20,423 22,199 23,399 25,271 26,891 28,551 30,615 32,039 34,199 36,479 37,631 38,807 41,579 46,619 50,171 51,527 52,895 55,215 57,119 59,999 63,999 67,071 70,215 72,359 74,519 77,279 78,959 82,343 89,351 94,859 96,719 98,591 104,279 110,879 116,255 120,407 122,495 126,015 131,027 136,151 140,615 144,395 148,215 153,647 158,399 163,199 170,543 175,559 180,599 185,759 189,215 193,595 198,015 204,287 209,759 212,519 215,291 222,747 232,307 238,139 244,019 249,995 255,015 264,159 271,439 281,879 294,839 303,575 312,471 319,215 323,759 328,319 337,535 346,911 354,015 358,799 363,599 370,871 376,991 380,687 389,339 403,199 410,879 414,731 421,191 429,015 434,279 443,519 454,271 461,031 470,579 482,999 495,599 508,343 521,267 531,431 540,215 547,595 556,499 566,999 574,559 583,679 592,895 606,791 625,655 643,167 654,479 664,199 674,039 678,971 683,927 693,863 713,975 729,311 734,447 739,595 755,111 770,879 776,159 781,451 802,715 824,459 835,379 851,903 868,607 879,839 889,239 900,591 919,631 937,019 946,719 958,431 972,179 986,039 Found 167 Frobenius numbers of consecutive primes under 1,000,000, of which 25 were under 10,000
REXX
<lang rexx>/*REXX program finds Frobenius numbers where the numbers are less than some number N. */ parse arg hi cols . /*obtain optional argument from the CL.*/ if hi== | hi=="," then hi= 10000 /* " " " " " " */ if cols== | cols=="," then cols= 10 /* " " " " " " */ call genP /*build array of semaphores for primes.*/ w= 10 /*width of a number in any column. */
@Frob= ' Frobenius numbers that are < ' commas(hi)
if cols>0 then say ' index │'center(@Frob, 1 + cols*(w+1) ) if cols>0 then say '───────┼'center("" , 1 + cols*(w+1), '─') idx= 1 /*initialize the index of output lines.*/ $= /*a list of Frobenius numbers (so far)*/
do j=1; jp= j + 1 /*generate Frobenius numbers < HI */
y= @.j * @.jp - @.j - @.jp
if y>= hi then leave
if cols==0 then iterate /*Build the list (to be shown later)? */
c= commas(y) /*maybe add commas to the number. */
$= $ right(c, max(w, length(c) ) ) /*add a Frobenius #──►list, allow big #*/
if j//cols\==0 then iterate /*have we populated a line of output? */
say center(idx, 7)'│' substr($, 2); $= /*display what we have so far (cols). */
idx= idx + cols /*bump the index count for the output*/
end /*j*/
if $\== then say center(idx, 7)"│" substr($, 2) /*possible display residual output.*/ say say 'Found ' commas(j-1) @FROB exit 0 /*stick a fork in it, we're all done. */ /*──────────────────────────────────────────────────────────────────────────────────────*/ commas: parse arg ?; do jc=length(?)-3 to 1 by -3; ?=insert(',', ?, jc); end; return ? /*──────────────────────────────────────────────────────────────────────────────────────*/ genP: @.1=2; @.2=3; @.3=5; @.4=7; @.5=11 /*define some low primes. */
#=5; s.#= @.# **2 /*number of primes so far; prime². */
/* [↓] generate more primes ≤ high.*/
do j=@.#+2 by 2 to hi /*find odd primes from here on. */
parse var j -1 _; if _==5 then iterate /*J divisible by 5? (right dig)*/
if j// 3==0 then iterate /*" " " 3? */
if j// 7==0 then iterate /*" " " 7? */
/* [↑] the above five lines saves time*/
do k=5 while s.k<=j /* [↓] divide by the known odd primes.*/
if j // @.k == 0 then iterate j /*Is J ÷ X? Then not prime. ___ */
end /*k*/ /* [↑] only process numbers ≤ √ J */
#= #+1; @.#= j; s.#= j*j /*bump # Ps; assign next P; P squared*/
end /*j*/; return</lang>
- output when using the default inputs:
index │ Frobenius numbers that are < 10,000 ───────┼─────────────────────────────────────────────────────────────────────────────────────────────────────────────── 1 │ 1 7 23 59 119 191 287 395 615 839 11 │ 1,079 1,439 1,679 1,931 2,391 3,015 3,479 3,959 4,619 5,039 21 │ 5,615 6,395 7,215 8,447 9,599 Found 25 Frobenius numbers that are < 10,000
Ring
<lang ring> load "stdlib.ring"
decimals(0) see "working..." + nl see "Frobenius numbers are:" + nl
row = 0 limit1 = 101 Frob = []
for n = 1 to limit1
if isprime(n)
add(Frob,n)
ok
next
for n = 1 to len(Frob)-1
row = row + 1
fr = Frob[n]*Frob[n+1]-Frob[n]-Frob[n+1]
see "" + fr + " "
if row%5 = 0
see nl
ok
next
see "Found " + row + " Frobenius primes" + nl see "done..." + nl </lang>
- Output:
working... Frobenius primes are: 1 7 23 59 119 191 287 395 615 839 1079 1439 1679 1931 2391 3015 3479 3959 4619 5039 5615 6395 7215 8447 9599 Found 25 Frobenius primes done...