Long primes: Difference between revisions
(promoted to full task status.) |
(Added Sidef) |
||
Line 933: | Line 933: | ||
return _</lang> |
return _</lang> |
||
{{out|output|text= is identical to the 1<sup>st</sup> REXX version.}} <br><br> |
{{out|output|text= is identical to the 1<sup>st</sup> REXX version.}} <br><br> |
||
=={{header|Sidef}}== |
|||
The smallest divisor d of p-1 such that 10^d = 1 (mod p), is the length of the period of the decimal expansion of 1/p. |
|||
<lang ruby>func is_long_prime(p) { |
|||
for d in (divisors(p-1)) { |
|||
if (powmod(10, d, p) == 1) { |
|||
return (d+1 == p) |
|||
} |
|||
} |
|||
return false |
|||
} |
|||
say "Long primes ≤ 500:" |
|||
say primes(500).grep(is_long_prime).join(' ') |
|||
for n in ([500, 1000, 2000, 4000, 8000, 16000, 32000, 64000]) { |
|||
say ("Number of long primes ≤ #{n}: ", primes(n).count_by(is_long_prime)) |
|||
}</lang> |
|||
{{out}} |
|||
<pre> |
|||
Long primes ≤ 500: |
|||
7 17 19 23 29 47 59 61 97 109 113 131 149 167 179 181 193 223 229 233 257 263 269 313 337 367 379 383 389 419 433 461 487 491 499 |
|||
Number of long primes ≤ 500: 35 |
|||
Number of long primes ≤ 1000: 60 |
|||
Number of long primes ≤ 2000: 116 |
|||
Number of long primes ≤ 4000: 218 |
|||
Number of long primes ≤ 8000: 390 |
|||
Number of long primes ≤ 16000: 716 |
|||
Number of long primes ≤ 32000: 1300 |
|||
Number of long primes ≤ 64000: 2430 |
|||
</pre> |
|||
=={{header|zkl}}== |
=={{header|zkl}}== |
Revision as of 08:05, 13 August 2018
You are encouraged to solve this task according to the task description, using any language you may know.
A long prime (the definition that will be used
here) are primes whose reciprocals (in decimal) have
a period length of one less than
the prime number (also expressed in decimal).
Long primes are also known as:
- base ten cyclic numbers
- full reptend primes
- golden primes
- long period primes
- maximal period primes
- proper primes
- Example
7 is the first long prime, the reciprocal of seven is 1/7, which is equal to the repeating decimal fraction 0.142857142857···
The length of the repeating part of the decimal fraction
is six, (the underlined part) which is one less
than the (decimal) prime number 7.
Thus 7 is a long prime.
There are other (more) general definitions of a long prime which
include wording/verbiage for other bases other than ten.
- Task
-
- Show all long primes up to 500 (preferably on one line).
- Show the number of long primes up to 500
- Show the number of long primes up to 1,000
- Show the number of long primes up to 2,000
- Show the number of long primes up to 4,000
- Show the number of long primes up to 8,000
- Show the number of long primes up to 16,000
- Show the number of long primes up to 32,000
- Show the number of long primes up to 64,000 (optional)
- Show all output here.
- Also see
- the Wikipedia webpage: full reptend prime.
- the MathWorld webpage: full reptend prime.
- the OEIS webpage: A1913.
C
<lang c>#include <stdio.h>
- include <stdlib.h>
- define TRUE 1
- define FALSE 0
typedef int bool;
void sieve(int limit, int primes[], int *count) {
bool *c = calloc(limit + 1, sizeof(bool)); // composite = TRUE // no need to process even numbers int i, p = 3, p2, n = 0; while (TRUE) { p2 = p * p; if (p2 > limit) break; for (i = p2; i <= limit; i += 2 * p) c[i] = TRUE; while (TRUE) { p += 2; if (!c[p]) break; } } for (i = 3; i <= limit; i += 2) { if (!c[i]) primes[n++] = i; } *count = n; free(c);
}
// finds the period of the reciprocal of n int findPeriod(int n) {
int i, r = 1, rr, period = 0; for (i = 1; i <= n + 1; ++i) { r = (10 * r) % n; } rr = r; while (TRUE) { r = (10 * r) % n; period++; if (r == rr) break; } return period;
}
int main() {
int i, prime, count = 0, index = 0, primeCount, longCount = 0; int *primes, *longPrimes; int numbers[] = {500, 1000, 2000, 4000, 8000, 16000, 32000, 64000}; int totals[8]; primes = calloc(6500, sizeof(int)); longPrimes = calloc(2500, sizeof(int)); sieve(64000, primes, &primeCount); for (i = 0; i < primeCount; ++i) { prime = primes[i]; if (findPeriod(prime) == prime - 1) { longPrimes[longCount++] = prime; } } for (i = 0; i < longCount; ++i) { if (longPrimes[i] > numbers[index]) { totals[index++] = count; } count++; } totals[7] = count; printf("The long primes up to 500 are:\n"); printf("["); for (i = 0; i < totals[0]; ++i) { printf("%d ", longPrimes[i]); } printf("\b]\n");
printf("\nThe number of long primes up to:\n"); for (i = 0; i < 8; ++i) { printf(" %5d is %d\n", numbers[i], totals[i]); } free(longPrimes); free(primes); return 0;
}</lang>
- Output:
The long primes up to 500 are: [7 17 19 23 29 47 59 61 97 109 113 131 149 167 179 181 193 223 229 233 257 263 269 313 337 367 379 383 389 419 433 461 487 491 499] The number of long primes up to: 500 is 35 1000 is 60 2000 is 116 4000 is 218 8000 is 390 16000 is 716 32000 is 1300 64000 is 2430
Go
<lang go>package main
import "fmt"
func sieve(limit int) []int {
var primes []int c := make([]bool, limit+1) // composite = true // no need to process even numbers p := 3 for { p2 := p * p if p2 > limit { break } for i := p2; i <= limit; i += 2 * p { c[i] = true } for { p += 2 if !c[p] { break } } } for i := 3; i <= limit; i += 2 { if !c[i] { primes = append(primes, i) } } return primes
}
// finds the period of the reciprocal of n func findPeriod(n int) int {
r := 1 for i := 1; i <= n+1; i++ { r = (10 * r) % n } rr := r period := 0 for { r = (10 * r) % n period++ if r == rr { break } } return period
}
func main() {
primes := sieve(64000) var longPrimes []int for _, prime := range primes { if findPeriod(prime) == prime-1 { longPrimes = append(longPrimes, prime) } } numbers := []int{500, 1000, 2000, 4000, 8000, 16000, 32000, 64000} index := 0 count := 0 totals := make([]int, len(numbers)) for _, longPrime := range longPrimes { if longPrime > numbers[index] { totals[index] = count index++ } count++ } totals[len(numbers)-1] = count fmt.Println("The long primes up to 500 are: ") fmt.Println(longPrimes[:totals[0]])
fmt.Println("\nThe number of long primes up to: ") for i, total := range totals { fmt.Printf(" %5d is %d\n", numbers[i], total) }
}</lang>
- Output:
The long primes up to 500 are: [7 17 19 23 29 47 59 61 97 109 113 131 149 167 179 181 193 223 229 233 257 263 269 313 337 367 379 383 389 419 433 461 487 491 499] The number of long primes up to: 500 is 35 1000 is 60 2000 is 116 4000 is 218 8000 is 390 16000 is 716 32000 is 1300 64000 is 2430
Kotlin
<lang scala>// Version 1.2.60
fun sieve(limit: Int): List<Int> {
val primes = mutableListOf<Int>() val c = BooleanArray(limit + 1) // composite = true // no need to process even numbers var p = 3 while (true) { val p2 = p * p if (p2 > limit) break for (i in p2..limit step 2 * p) c[i] = true while (true) { p += 2 if (!c[p]) break } } for (i in 3..limit step 2) { if (!c[i]) primes.add(i) } return primes
}
// finds the period of the reciprocal of n fun findPeriod(n: Int): Int {
var r = 1 for (i in 1..n + 1) r = (10 * r) % n val rr = r var period = 0 while (true) { r = (10 * r) % n period++ if (r == rr) break } return period
}
fun main(args: Array<String>) {
val primes = sieve(64000) val longPrimes = mutableListOf<Int>() for (prime in primes) { if (findPeriod(prime) == prime - 1) { longPrimes.add(prime) } } val numbers = listOf(500, 1000, 2000, 4000, 8000, 16000, 32000, 64000) var index = 0 var count = 0 val totals = IntArray(numbers.size) for (longPrime in longPrimes) { if (longPrime > numbers[index]) { totals[index++] = count } count++ } totals[numbers.lastIndex] = count println("The long primes up to 500 are:") println(longPrimes.take(totals[0]))
println("\nThe number of long primes up to:") for ((i, total) in totals.withIndex()) { System.out.printf(" %5d is %d\n", numbers[i], total) }
}</lang>
- Output:
The long primes up to 500 are: [7, 17, 19, 23, 29, 47, 59, 61, 97, 109, 113, 131, 149, 167, 179, 181, 193, 223, 229, 233, 257, 263, 269, 313, 337, 367, 379, 383, 389, 419, 433, 461, 487, 491, 499] The number of long primes up to: 500 is 35 1000 is 60 2000 is 116 4000 is 218 8000 is 390 16000 is 716 32000 is 1300 64000 is 2430
Pascal
first post.old program modified. Using Euler Phi
www . arndt-bruenner.de/mathe/scripts/periodenlaenge.htm
<lang pascal> PROGRAM Periode;
{$IFDEF FPC}
{$MODE Delphi} {$OPTIMIZATION ON} {$OPTIMIZATION Regvar} {$OPTIMIZATION Peephole} {$OPTIMIZATION cse} {$OPTIMIZATION asmcse}
{$else}
{$Apptype Console}
{$ENDIF}
uses
sysutils;
const
cBASIS = 10; PRIMFELDOBERGRENZE = 6542; {Das sind alle Primzahlen bis 2^16} {Das reicht fuer al8le Primzahlen bis 2^32} TESTZAHL = 500;//429496709;//High(Dword) DIV cBasis;
type
tPrimFeld = array [1..PRIMFELDOBERGRENZE] of Word; tFaktorPotenz = record Faktor, Potenz : DWord; end; //2*3*5*7*11*13*17*19*23 *29 > DWord also maximal 9 Faktoren tFaktorFeld = array [1..9] of TFaktorPotenz;//DWord
// tFaktorFeld = array [1..15] of TFaktorPotenz;//QWord
tFaktorisieren = class(TObject) private fFakZahl : DWord; fFakBasis : DWord; fFakAnzahl : Dword; fAnzahlMoeglicherTeiler : Dword; fEulerPhi : DWORD; fStartPeriode : DWORD; fPeriodenLaenge : DWORD; fTeiler : array of DWord; fFaktoren : tFaktorFeld; fBasFakt : tFaktorFeld; fPrimfeld : tPrimFeld;
procedure PrimFeldAufbauen; procedure Fakteinfuegen(var Zahl:Dword;inFak:Dword); function BasisPeriodeExtrahieren(var inZahl:Dword):DWord; procedure NachkommaPeriode(var OutText: String); public constructor create; overload; function Prim(inZahl:Dword):Boolean; procedure AusgabeFaktorfeld(n : DWord); procedure Faktorisierung (inZahl: DWord); procedure TeilerErmitteln; procedure PeriodeErmitteln(inZahl:Dword); function BasExpMod( b, e, m : Dword) : DWord; property EulerPhi : Dword read fEulerPhi; property PeriodenLaenge: DWord read fPeriodenLaenge ; property StartPeriode: DWord read fStartPeriode ; end;
constructor tFaktorisieren.create; begin
inherited; PrimFeldAufbauen; fFakZahl := 0; fFakBasis := cBASIS; Faktorisierung(fFakBasis); fBasFakt := fFaktoren;
fFakZahl := 0; fEulerPhi := 1; fPeriodenLaenge :=0; fFakZahl := 0; fFakAnzahl := 0; fAnzahlMoeglicherTeiler := 0;
end;
function tFaktorisieren.Prim(inZahl:Dword):Boolean; {Testet auf PrimZahl} var
Wurzel, pos : Dword;
Begin
if fFakZahl = inZahl then begin result := (fAnzahlMoeglicherTeiler = 2); exit; end; result := false; if inZahl >1 then begin result := true; Pos := 1; Wurzel:= trunc(sqrt(inZahl)); While fPrimFeld[Pos] <= Wurzel do begin if (inZahl mod fPrimFeld[Pos])=0 then begin result := false; break; end; inc(Pos); IF Pos > High(fPrimFeld) then break; end; end;
end;
Procedure tFaktorisieren.PrimFeldAufbauen; {Baut die Liste der Primzahlen bis Obergrenze auf} const
MAX = 65536;
var
TestaufPrim, Zaehler,delta : Dword;
begin Zaehler := 1; fPrimFeld[Zaehler] := 2; inc(Zaehler); fPrimFeld[Zaehler] := 3;
delta := 2; TestAufPrim:=5; repeat
if prim(TestAufPrim) then begin inc(Zaehler); fPrimFeld[Zaehler] := TestAufPrim; end; inc(TestAufPrim,delta); delta := 6-delta; // 2,4,2,4,2,4,2,
until (TestAufPrim>=MAX);
end; {PrimfeldAufbauen}
procedure tFaktorisieren.Fakteinfuegen(var Zahl:Dword;inFak:Dword);
var
i : DWord;
begin
inc(fFakAnzahl); with fFaktoren[fFakAnzahl] do begin fEulerPhi := fEulerPhi*(inFak-1); Faktor :=inFak; Potenz := 0; while (Zahl mod inFak) = 0 do begin Zahl := Zahl div inFak; inc(Potenz); end; For i := 2 to Potenz do fEulerPhi := fEulerPhi*inFak; end; fAnzahlMoeglicherTeiler:=fAnzahlMoeglicherTeiler*(1+fFaktoren[fFakAnzahl].Potenz);
end;
procedure tFaktorisieren.Faktorisierung (inZahl: DWord); var
j, og : longint;
begin if fFakZahl = inZahl then
exit;
fPeriodenLaenge := 0; fFakZahl := inZahl; fEulerPhi := 1; fFakAnzahl := 0; fAnzahlMoeglicherTeiler :=1; setlength(fTeiler,0);
If inZahl < 2 then
exit;
og := round(sqrt(inZahl)+1.0); {Suche Teiler von inZahl} for j := 1 to High(fPrimfeld) do
begin If fPrimfeld[j]> OG then Break; if (inZahl mod fPrimfeld[j])= 0 then Fakteinfuegen(inZahl,fPrimfeld[j]); end;
If inZahl>1 then
Fakteinfuegen(inZahl,inZahl);
TeilerErmitteln; end; {Faktorisierung}
procedure tFaktorisieren.AusgabeFaktorfeld(n : DWord); var
i :integer;
begin
if fFakZahl <> n then Faktorisierung(n); write(fAnzahlMoeglicherTeiler:5,' Faktoren ');
For i := 1 to fFakAnzahl-1 do with fFaktoren[i] do IF potenz >1 then write(Faktor,'^',Potenz,'*') else write(Faktor,'*'); with fFaktoren[fFakAnzahl] do IF potenz >1 then write(Faktor,'^',Potenz) else write(Faktor);
writeln(' Euler Phi: ',fEulerPhi:12,PeriodenLaenge:12);
end;
procedure tFaktorisieren.TeilerErmitteln; var
Position : DWord; i,j : DWord; procedure FaktorAufbauen(Faktor: DWord;n: DWord); var i, Pot : DWord; begin Pot := 1; i := 0; repeat IF n > Low(fFaktoren) then FaktorAufbauen(Pot*Faktor,n-1) else begin FTeiler[Position] := Pot*Faktor; inc(Position); end; Pot := Pot*fFaktoren[n].Faktor; inc(i); until i > fFaktoren[n].Potenz; end;
begin
Position:= 0; setlength(FTeiler,fAnzahlMoeglicherTeiler); FaktorAufbauen(1,fFakAnzahl); //Sortieren For i := Low(fTeiler) to fAnzahlMoeglicherTeiler-2 do begin j := i; while (j>=Low(fTeiler)) AND (fTeiler[j]>fTeiler[j+1]) do begin Position := fTeiler[j]; fTeiler[j] := fTeiler[j+1]; fTeiler[j+1]:= Position; dec(j); end; end;
end;
function tFaktorisieren.BasisPeriodeExtrahieren(var inZahl:Dword):DWord; var
i,cnt, Teiler: Dword;
begin
cnt := 0; result := 0; For i := Low(fBasFakt) to High(fBasFakt) do begin with fBasFakt[i] do begin IF Faktor = 0 then BREAK; Teiler := Faktor; For cnt := 2 to Potenz do Teiler := Teiler*Faktor; end; cnt := 0; while (inZahl<> 0) AND (inZahl mod Teiler = 0) do begin inZahl := inZahl div Teiler; inc(cnt); end; IF cnt > result then result := cnt; end;
end;
procedure tFaktorisieren.PeriodeErmitteln(inZahl:Dword); var
i, TempZahl, TempPhi, TempPer, TempBasPer: DWord;
begin
Faktorisierung(inZahl); TempZahl := inZahl; //Die Basis_Nicht_Periode ermitteln TempBasPer := BasisPeriodeExtrahieren(TempZahl); TempPer := 0; IF TempZahl >1 then begin Faktorisierung(TempZahl); TempPhi := fEulerPhi; IF (TempPhi > 1) then begin Faktorisierung(TempPhi); i := 0; repeat TempPer := fTeiler[i]; IF BasExpMod(fFakBasis,TempPer,TempZahl)= 1 then Break; inc(i); until i >= Length(fTeiler); IF i >= Length(fTeiler) then TempPer := inZahl-1; end; end;
Faktorisierung(inZahl); fPeriodenlaenge := TempPer; fStartPeriode := TempBasPer;
end;
procedure tFaktorisieren.NachkommaPeriode(var OutText: String); var
i, limit : integer; Rest, Rest1, Divi, basis: DWord; pText : pChar; procedure Ziffernfolge(Ende: longint); var j : longint; begin j := i-Ende; while j < 0 do begin Rest := Rest *Basis; Rest1:= Rest Div Divi; Rest := Rest-Rest1*Divi;//== Rest1 Mod Divi pText^ := chr(Rest1+Ord('0')); inc(pText);
inc(j);
end;
i := Ende;
end;
begin
limit:= fStartPeriode+fPeriodenlaenge; setlength(OutText,limit+2+2+5); OutText[1]:='0'; OutText[2]:='.'; pText := @OutText[3];
Rest := 1; Divi := fFakZahl; Basis := fFakBasis; i := 0; Ziffernfolge(fStartPeriode); if fPeriodenlaenge = 0 then begin setlength(OutText,fStartPeriode+2); EXIT; end;
pText^ := '_'; inc(pText); Ziffernfolge(limit); pText^ := '_'; inc(pText);
Ziffernfolge(limit+5);
end;
type
tZahl = integer; tRestFeld = array[0..31] of integer;
VAR
F : tFaktorisieren;
function tFaktorisieren.BasExpMod( b, e, m : Dword) : DWord; begin
Result := 1; IF m = 0 then exit; Result := 1; while ( e > 0 ) do begin if (e AND 1) <> 0 then Result := (Result * int64(b)) mod m; b := (int64(b) * b ) mod m; e := e shr 1; end;
end;
procedure start; VAR
Limit, Testzahl : DWord; longPrimCount : int64; t1,t0: TDateTime;
BEGIN
Limit := 500; Testzahl := 2; longPrimCount := 0; t0 := time;
repeat write(Limit:8,': '); repeat if F.Prim(Testzahl) then begin F.PeriodeErmitteln(Testzahl); if F.PeriodenLaenge = Testzahl-1 then Begin inc(longPrimCount); IF Limit = 500 then write(TestZahl,','); end end; inc(Testzahl); until TestZahl = Limit; inc(Limit,Limit); write(' .. count ',longPrimCount:8,' '); t1:= time; If (t1-t0)>1/864000 then write(FormatDateTime('HH:NN:SS.ZZZ',t1-T0)); writeln; until Limit > 10*1000*1000;
t1 := time; writeln; writeln('count of long primes ',longPrimCount); writeln('Benoetigte Zeit ',FormatDateTime('HH:NN:SS.ZZZ',T1-T0));
END;
BEGIN
F := tFaktorisieren.create; writeln('Start'); start; writeln('Fertig.'); F.free; readln;
end.</lang>
- Output:
sh-4.4# ./Periode Start 500: 7,17,19,23,29,47,59,61,97,109,113,131,149,167,179,181,193,223,229,233,257,263,269,313,337,367,379,383,389,419,433,461,487,491,499, .. count 35 1000: .. count 60 2000: .. count 116 4000: .. count 218 8000: .. count 390 16000: .. count 716 32000: .. count 1300 64000: .. count 2430 128000: .. count 4498 256000: .. count 8434 00:00:00.100 512000: .. count 15920 00:00:00.220 1024000: .. count 30171 00:00:00.494 2048000: .. count 57115 00:00:01.140 4096000: .. count 108381 00:00:02.578 8192000: .. count 206594 00:00:06.073 count of long primes 206594 Benoetigte Zeit 00:00:06.073 Fertig.
Perl 6
Not very fast as the numbers get larger. 64000 takes a little over 15 minutes on my computer. 😕 <lang perl6>my @long-primes = lazy (1..*).grep(*.is-prime).hyper(:8degree, :8batch).grep({1+(1/$_).base-repeating[1].chars == $_});
put "Long primes ≤ 500:\n", @long-primes[^(@long-primes.first: * > 500, :k)];
say "\nNumber of long primes ≤ $_: ", +@long-primes[^(@long-primes.first: * > $_, :k)]
for 500, 1000, 2000, 4000, 8000, 16000, 32000, 64000;</lang>
- Output:
Long primes ≤ 500: 7 17 19 23 29 47 59 61 97 109 113 131 149 167 179 181 193 223 229 233 257 263 269 313 337 367 379 383 389 419 433 461 487 491 499 Number of long primes ≤ 500: 35 Number of long primes ≤ 1000: 60 Number of long primes ≤ 2000: 116 Number of long primes ≤ 4000: 218 Number of long primes ≤ 8000: 390 Number of long primes ≤ 16000: 716 Number of long primes ≤ 32000: 1300 Number of long primes ≤ 64000: 2430
REXX
For every doubling of the limit, it takes about roughly 8 times longer to compute the long primes.
uses odd numbers
<lang rexx>/*REXX pgm calculates/displays base ten long primes (AKA golden primes, proper primes,*/ /*───────────────────── maximal period primes, long period primes, full reptend primes).*/ parse arg a /*obtain optional argument from the CL.*/ if a= | a="," then a= '500 -500 -100 -2000 -4000 -8000 -16000' /*use the default?*/
do k=1 for words(a); H=word(a, k) /*step through the list of high limits.*/ neg= H<1 /*used as an indicator to display count*/ H= abs(H) /*obtain the absolute value of H. */ numeric digits max(H, 500) /*insure enough dec digs for periodLen.*/ $= /*the list of long primes (so far). */ do j=7 to H by 2 /*start with 7, just use odd integers.*/ if .len(j) + 1 \== j then iterate /*period length too small? " " " */ $=$ j /*add the long prime to the $ list.*/ end /*j*/ say if neg then do; say 'number of long primes ≤ ' H " is: " words($); end else do; say 'list of long primes ≤ ' H":"; say strip($); end end /*k*/
exit /*stick a fork in it, we're all done. */ /*──────────────────────────────────────────────────────────────────────────────────────*/ .len: procedure; parse arg x -1 z; y=9 /*obtain the argument from the caller. */
if z==5 then return 0 /*if the last digit is 5, then skip. */ _=1 do while y//x \== 0; y= y'9'; _= length(y) end /*while*/ return _</lang>
list of long primes ≤ 500: 7 17 19 23 29 47 59 61 97 109 113 131 149 167 179 181 193 223 229 233 257 263 269 313 337 367 379 383 389 419 433 461 487 491 499 number of long primes ≤ 500 is: 35 number of long primes ≤ 1000 is: 60 number of long primes ≤ 2000 is: 116 number of long primes ≤ 4000 is: 218 number of long primes ≤ 8000 is: 390 number of long primes ≤ 16000 is: 716 number of long primes ≤ 32000 is: 1300
uses primes
This REXX version is about 15% faster than the 1st REXX version (becauses it only tests primes). <lang rexx>/*REXX pgm calculates/displays base ten long primes (AKA golden primes, proper primes,*/ /*───────────────────── maximal period primes, long period primes, full reptend primes).*/ parse arg a /*obtain optional argument from the CL.*/ if a= | a="," then a= '500 -500 -1000 -2000 -4000 -8000 -16000 -32000' /*use default?*/ m=0; aa=words(a) /* [↑] two list types of low primes. */
do j=1 for aa; m= max(m, abs(word(a, j))) /*find the maximum argument in the list*/ end /*j*/
call genP /*go and generate some primes. */
do k=1 for aa; H=word(a, k) /*step through the list of high limits.*/ neg= H<1 /*used as an indicator to display count*/ H= abs(H) /*obtain the absolute value of H. */ numeric digits max(H, 500) /*insure enough dec digs for periodLen.*/ $= /*the list of long primes (so far). */ do j=7 to H by 2 if \@.j then iterate /*Is J not a prime? Then skip it. */ if .len(j) + 1 \== j then iterate /*period length too small? " " " */ $=$ j /*add the long prime to the $ list.*/ end /*j*/ /* [↑] some pretty weak prime testing.*/ say if neg then do; say 'number of long primes ≤ ' H " is: " words($); end else do; say 'list of long primes ≤ ' H":"; say strip($); end end /*k*/
exit /*stick a fork in it, we're all done. */ /*──────────────────────────────────────────────────────────────────────────────────────*/ genP: @.=0; @.2=1; @.3=1; @.5=1; @.7=1; @.11=1; !.=0; !.1=2; !.2=3; !.3=5; !.4=7; !.5=11
#=5 /*the number of primes (so far). */ do g=!.#+2 by 2 until #>=m /*gen enough primes to satisfy max A. */ do d=2 until !.d**2 > g /*only divide up to square root of X. */ if g // !.d == 0 then iterate g /*Divisible? Then skip this integer. */ end /*d*/ /* [↓] a spanking new prime was found.*/ #=#+1; @.g=1; !.#=g /*bump P counter; assign P, add to P's.*/ end /*g*/ return
/*──────────────────────────────────────────────────────────────────────────────────────*/ .len: procedure; parse arg x; _=1; y=9 /*obtain the argument from the caller. */
do while y//x \== 0; y= y'9'; _= length(y) end /*while*/ return _</lang>
- output is identical to the 1st REXX version.
Sidef
The smallest divisor d of p-1 such that 10^d = 1 (mod p), is the length of the period of the decimal expansion of 1/p. <lang ruby>func is_long_prime(p) {
for d in (divisors(p-1)) { if (powmod(10, d, p) == 1) { return (d+1 == p) } }
return false
}
say "Long primes ≤ 500:" say primes(500).grep(is_long_prime).join(' ')
for n in ([500, 1000, 2000, 4000, 8000, 16000, 32000, 64000]) {
say ("Number of long primes ≤ #{n}: ", primes(n).count_by(is_long_prime))
}</lang>
- Output:
Long primes ≤ 500: 7 17 19 23 29 47 59 61 97 109 113 131 149 167 179 181 193 223 229 233 257 263 269 313 337 367 379 383 389 419 433 461 487 491 499 Number of long primes ≤ 500: 35 Number of long primes ≤ 1000: 60 Number of long primes ≤ 2000: 116 Number of long primes ≤ 4000: 218 Number of long primes ≤ 8000: 390 Number of long primes ≤ 16000: 716 Number of long primes ≤ 32000: 1300 Number of long primes ≤ 64000: 2430
zkl
Using GMP (GNU Multiple Precision Arithmetic Library, probabilistic primes), because it is easy and fast to generate primes. <lang zkl>var [const] BN=Import("zklBigNum"); // libGMP primes,p := List.createLong(7_000), BN(3); // one big alloc vs lots of allocs while(p.nextPrime()<=64_000){ primes.append(p.toInt()) } // 6412 of them, skipped 2 primes.append(p.toInt()); // and one more so tail prime is >64_000
fcn findPeriod(n){
r,period := 1,0; do(n){ r=(10*r)%n } rr:=r; while(True){ // reduce is more concise but 2.5 times slower r=(10*r)%n; period+=1; if(r==rr) break; } period
}</lang> <lang zkl>longPrimes:=primes.filter(fcn(p){ findPeriod(p)==p-1 }); // yawn fiveHundred:=longPrimes.filter('<(500)); println("The long primes up to 500 are:\n",longPrimes.filter('<(500)).concat(","));
println("\nThe number of long primes up to:"); foreach n in (T(500, 1000, 2000, 4000, 8000, 16000, 32000, 64000)){
println(" %5d is %d".fmt( n, longPrimes.filter1n('>(n)) ));
}</lang>
- Output:
The long primes up to 500 are: 7,17,19,23,29,47,59,61,97,109,113,131,149,167,179,181,193,223,229,233,257,263,269,313,337,367,379,383,389,419,433,461,487,491,499 The number of long primes up to: 500 is 35 1000 is 60 2000 is 116 4000 is 218 8000 is 390 16000 is 716 32000 is 1300 64000 is 2430