Wagstaff primes: Difference between revisions
m
→{{header|Wren}}: Changed to Wren S/H
(→{{header|J}}: added Java version) |
m (→{{header|Wren}}: Changed to Wren S/H) |
||
| (33 intermediate revisions by 18 users not shown) | |||
Line 1:
{{
;Definition
A ''Wagstaff prime'' is a prime number of the form ''(2^p + 1)/3'' where the exponent ''p'' is an odd prime.
Line 83:
</pre>
=={{header|ALGOL W}}==
As Algol W integers are limited to 32 bits, 64 bit reals are used which can accurately represent integers up to 2^53.
<syntaxhighlight lang="algolw">
begin % find some Wagstaff primes: primes of the form ( 2^p + 1 ) / 3 %
% where p is an odd prime %
integer wCount;
% returns true if d exactly divides v, false otherwise %
logical procedure divides( long real value d, v ) ;
begin
long real q, p10;
q := v / d;
p10 := 1;
while p10 * 10 < q do begin
p10 := p10 * 10
end while_p10_lt_q ;
while p10 >= 1 do begin
while q >= p10 do q := q - p10;
p10 := p10 / 10
end while_p10_ge_1 ;
q = 0
end divides ;
% prints v as a 14 digit integer number %
procedure printI14( long real value v ) ;
begin
integer iv;
long real r;
r := abs( v );
if v < 0 then writeon( s_w := 0, "-" );
iv := truncate( r / 1000000 );
if iv < 1 then begin
writeon( i_w := 8, s_w := 0, " ", truncate( r ) )
end
else begin
string(6) sValue;
writeon( i_w := 8, s_w := 0, iv );
iv := truncate( abs( r ) - ( iv * 1000000.0 ) );
for sPos := 5 step -1 until 0 do begin
sValue( sPos // 1 ) := code( ( iv rem 10 ) + decode( "0" ) );
iv := iv div 10
end for_sPos;
writeon( s_w := 0, sValue )
end if_iv_lt_1__
end printI ;
wCount := 0;
% find the Wagstaff primes using long reals to hold the numbers, which %
% accurately represent integers up to 2^53, so only consider primes < 53 %
for p := 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47 do begin
long real w, powerOfTwo;
logical isPrime;
powerOfTwo := 1;
for i := 1 until p do powerOfTwo := powerOfTwo * 2;
w := ( powerOfTwo + 1 ) / 3;
isPrime := not divides( 2, w );
if isPrime then begin
integer f, toNext;
long real f2;
f := 3;
f2 := 9;
toNext := 16;
while isPrime and f2 <= w do begin
isPrime := not divides( f, w );
f := f + 2;
f2 := f2 + toNext;
toNext := toNext + 8
end while_isPrime_and_f2_le_x
end if_isPrime ;
if isPrime then begin
wCount := wCount + 1;
write( i_w := 3, s_w := 0, wCount, ": ", p, " " );
if p >= 32 then printI14( w )
else begin
integer iw;
iw := truncate( w );
writeon( i_w := 14, s_w := 0, iw )
end of_p_ge_32__ ;
if wCount >= 10 then goto done % stop at 10 Wagstaff primes %
end if_isPrime
end for_p ;
done:
end.
</syntaxhighlight>
{{out}}
<pre>
1: 3 3
2: 5 11
3: 7 43
4: 11 683
5: 13 2731
6: 17 43691
7: 19 174763
8: 23 2796203
9: 31 715827883
10: 43 2932031007403
</pre>
=={{header|Arturo}}==
<syntaxhighlight lang="arturo">wagstaff?: function [e][
and? -> prime? e -> prime? (1+2^e)/3
]
summarize: function [n][
n: ~"|n|"
s: size n
if s > 20 -> n: ((take n 10)++"...")++drop.times:s-10 n
n ++ ~" (|s| digits)"
]
exponents: select.first:24 range.step:2 1 ∞ => wagstaff?
loop.with:'i exponents 'x -> print [
pad ~"|i+1|:" 3 pad ~"|x| -" 6 summarize (1+2^x)/3
]</syntaxhighlight>
{{out}}
<pre>
p = 3, m = 3
p = 5, m = 11
p = 7, m = 43
p = 11, m = 683
p = 13, m = 2731
p = 17, m = 43691
p = 19, m = 174763
p = 23, m = 2796203
p = 31, m = 715827883
p = 43, m = 2932031007403
p = 61, m = 768614336404564651
p = 79, m = 2014876366...8195784363 (24 digits)
p = 101, m = 8451004001...1135470251 (30 digits)
p = 127, m = 5671372782...8628035243 (38 digits)
p = 167, m = 6235740319...3121833643 (50 digits)
p = 191, m = 1046183622...7339085483 (58 digits)
p = 199, m = 2678230073...8805883563 (60 digits)
p = 313, m = 5562466239...8130434731 (94 digits)
p = 347, m = 9556244233...1903606443 (104 digits)
p = 701, m = 3506757267...7823854251 (211 digits)
p = 1709, m = 9619252724...9070528171 (514 digits)
p = 2617, m = 2081504709...3435947691 (788 digits)
p = 3539, m = 7379609820...6486497963 (1065 digits)
p = 5807, m = 4018496237...6663568043 (1748 digits)
</pre>
=={{header|BASIC256}}==
{{trans|FreeBASIC}}
<syntaxhighlight lang="BASIC256">function isPrime(v)
if v < 2 then return False
if v mod 2 = 0 then return v = 2
if v mod 3 = 0 then return v = 3
d = 5
while d * d <= v
if v mod d = 0 then return False else d += 2
end while
return True
end function
subroutine Wagstaff(num)
pri = 1
wcount = 0
wag = 0
while wcount < num
pri += 2
if isPrime(pri) then
wag = (2 ^ pri + 1) / 3
if isPrime(wag) then
wcount += 1
print rjust(wcount,2); ": "; rjust(pri,2); " => "; int(wag)
end if
end if
end while
end subroutine
call Wagstaff(9) #BASIC-256 does not allow larger numbers
end</syntaxhighlight>
=={{header|C}}==
{{libheader|GMP}}
This takes about 24 seconds to find the first 24 Wagstaff primes but 589 seconds to find the first 29.
<syntaxhighlight lang="c">#include <stdio.h>
#include <string.h>
#include <gmp.h>
int main() {
const int limit = 29;
int count = 0;
char tmp[40];
mpz_t p, w;
mpz_init_set_ui(p, 1);
mpz_init(w);
while (count < limit) {
mpz_nextprime(p, p);
mpz_set_ui(w, 1);
unsigned long ulp = mpz_get_ui(p);
mpz_mul_2exp(w, w, ulp);
mpz_add_ui(w, w, 1);
mpz_tdiv_q_ui(w, w, 3);
if (mpz_probab_prime_p(w, 15) > 0) {
++count;
char *ws = mpz_get_str(NULL, 10, w);
size_t le = strlen(ws);
if (le < 34) {
strcpy(tmp, ws);
} else {
strncpy(tmp, ws, 15);
strcpy(tmp + 15, "...");
strncpy(tmp + 18, ws + le - 15, 16);
}
printf("%5lu: %s", ulp, tmp);
if (le >=34) printf( " (%ld digits)", le);
printf("\n");
}
}
return 0;
}</syntaxhighlight>
{{out}}
<pre>
3: 3
5: 11
7: 43
11: 683
13: 2731
17: 43691
19: 174763
23: 2796203
31: 715827883
43: 2932031007403
61: 768614336404564651
79: 201487636602438195784363
101: 845100400152152934331135470251
127: 567137278201564...101238628035243 (38 digits)
167: 623574031927851...838653121833643 (50 digits)
191: 104618362256444...574077339085483 (58 digits)
199: 267823007376498...963798805883563 (60 digits)
313: 556246623937737...099018130434731 (94 digits)
347: 955624423329196...712921903606443 (104 digits)
701: 350675726769891...862167823854251 (211 digits)
1709: 961925272487010...857299070528171 (514 digits)
2617: 208150470990263...847933435947691 (788 digits)
3539: 737960982013072...304686486497963 (1065 digits)
5807: 401849623734300...466686663568043 (1748 digits)
10501: 435374724597165...340431558126251 (3161 digits)
10691: 683222859828090...209259644365483 (3218 digits)
11279: 692149388152358...081084271176363 (3395 digits)
12391: 385083595083686...967789680077483 (3730 digits)
14479: 136831460986221...678801865288363 (4359 digits)
</pre>
=={{header|C++}}==
{{libheader|GMP}}
{{libheader|Primesieve}}
<syntaxhighlight lang="cpp">#include <gmpxx.h>
#include <primesieve.hpp>
#include <iostream>
using big_int = mpz_class;
std::string to_string(const big_int& num, size_t max_digits) {
std::string str = num.get_str();
size_t len = str.size();
if (len > max_digits) {
str.replace(max_digits / 2, len - max_digits, "...");
str += " (";
str += std::to_string(len);
str += " digits)";
}
return str;
}
bool is_probably_prime(const big_int& n) {
return mpz_probab_prime_p(n.get_mpz_t(), 25) != 0;
}
int main() {
const big_int one(1);
primesieve::iterator pi;
pi.next_prime();
for (int i = 0; i < 24;) {
uint64_t p = pi.next_prime();
big_int n = ((one << p) + 1) / 3;
if (is_probably_prime(n))
std::cout << ++i << ": " << p << " - " << to_string(n, 30) << '\n';
}
}</syntaxhighlight>
{{out}}
<pre>
1: 3 - 3
2: 5 - 11
3: 7 - 43
4: 11 - 683
5: 13 - 2731
6: 17 - 43691
7: 19 - 174763
8: 23 - 2796203
9: 31 - 715827883
10: 43 - 2932031007403
11: 61 - 768614336404564651
12: 79 - 201487636602438195784363
13: 101 - 845100400152152934331135470251
14: 127 - 567137278201564...101238628035243 (38 digits)
15: 167 - 623574031927851...838653121833643 (50 digits)
16: 191 - 104618362256444...574077339085483 (58 digits)
17: 199 - 267823007376498...963798805883563 (60 digits)
18: 313 - 556246623937737...099018130434731 (94 digits)
19: 347 - 955624423329196...712921903606443 (104 digits)
20: 701 - 350675726769891...862167823854251 (211 digits)
21: 1709 - 961925272487010...857299070528171 (514 digits)
22: 2617 - 208150470990263...847933435947691 (788 digits)
23: 3539 - 737960982013072...304686486497963 (1065 digits)
24: 5807 - 401849623734300...466686663568043 (1748 digits)
</pre>
=={{header|Craft Basic}}==
<syntaxhighlight lang="basic">'Craft Basic may only calculate up to 9
let n = 9
let p = 1
do
let p = p + 2
if prime(p) then
let w = (2 ^ p + 1) / 3
if prime(w) then
let c = c + 1
print tab, c, tab, p, tab, w
endif
endif
title p
wait
loop c < n</syntaxhighlight>
{{out}}
<pre>
1 3 3
2 5 11
3 7 43
4 11 683
5 13 2731
6 17 43691
7 19 174763
8 23 2796203
9 31 715827883
</pre>
=={{header|Delphi}}==
{{works with|Delphi|6.0}}
{{libheader|SysUtils,StdCtrls}}
Finds Wagstaff primes up to 2^64, the biggest integer 32-bit Delphi supports.
<syntaxhighlight lang="Delphi">
procedure WagstaffPrimes(Memo: TMemo);
{Finds Wagstaff primes up to 6^64}
var P,B,R: int64;
var Cnt: integer;
begin
{B is int64 to force 64-bit arithmetic}
P:=0; B:=1; Cnt:=0;
Memo.Lines.Add(' #: P (2^P + 1)/3');
while P<64 do
begin
R:=((B shl P) + 1) div 3;
if IsPrime(P) and IsPrime(R) then
begin
Inc(Cnt);
Memo.Lines.Add(Format('%2d: %2d %24.0n',[Cnt,P,R+0.0]));
end;
Inc(P);
end;
end;
</syntaxhighlight>
{{out}}
<pre>
#: P (2^P + 1)/3
1: 3 3
2: 5 11
3: 7 43
4: 11 683
5: 13 2,731
6: 17 43,691
7: 19 174,763
8: 23 2,796,203
9: 31 715,827,883
10: 43 2,932,031,007,403
11: 61 768,614,336,404,564,651
Elapsed Time: 29.504 Sec.
</pre>
=={{header|EasyLang}}==
<syntaxhighlight lang="easylang">
func prime n .
if n mod 2 = 0 and n > 2
return 0
.
i = 3
while i <= sqrt n
if n mod i = 0
return 0
.
i += 2
.
return 1
.
pri = 1
while nwag <> 10
pri += 2
if prime pri = 1
wag = (pow 2 pri + 1) / 3
if prime wag = 1
nwag += 1
print pri & " => " & wag
.
.
.
</syntaxhighlight>
=={{header|F_Sharp|F#}}==
Line 107 ⟶ 545:
61 79 101 127 167 191 199 313 347 701 1709 2617 3539 5807
</pre>
=={{header|FreeBASIC}}==
{{trans|Python}}
<syntaxhighlight lang="freebasic">Function isPrime(Byval num As Ulongint) As Boolean
If num < 2 Then Return False
If num Mod 2 = 0 Then Return num = 2
If num Mod 3 = 0 Then Return num = 3
Dim d As Uinteger = 5
While d * d <= num
If num Mod d = 0 Then Return False Else d += 2
Wend
Return True
End Function
Sub Wagstaff(num As Ulongint)
Dim As Ulongint pri = 1, wcount = 0, wag
While wcount < num
pri += 2
If isPrime(pri) Then
wag = (2 ^ pri + 1) / 3
If isPrime(wag) Then
wcount += 1
Print Using "###: ### => ##,###,###,###,###"; wcount; pri; wag
End If
End If
Wend
End Sub
Wagstaff(10)
Sleep</syntaxhighlight>
{{out}}
<pre> 1: 3 => 3
2: 5 => 11
3: 7 => 43
4: 11 => 683
5: 13 => 2,731
6: 17 => 43,691
7: 19 => 174,763
8: 23 => 2,796,203
9: 31 => 715,827,883
10: 43 => 2,932,031,007,403</pre>
=={{header|Gambas}}==
{{trans|FreeBASIC}}
<syntaxhighlight lang="vbnet">Use "isprime.bas"
Public Sub Main()
Wagstaff(10)
End
Sub Wagstaff(num As Long)
Dim pri As Long = 1
Dim wcount As Long = 0
Dim wag As Long
While wcount < num
pri = pri + 2
If isPrime(pri) Then
wag = (2 ^ pri + 1) / 3
If isPrime(wag) Then
wcount += 1
Print Format$(Str(wcount), "###"); ": "; Format$(Str(pri), "###"); " => "; Int(wag)
End If
End If
Wend
End Sub</syntaxhighlight>
{{out}}
<pre>Similar to FreeBASIC entry.</pre>
=={{header|Go}}==
{{trans|C}}
{{libheader|GMP(Go wrapper)}}
{{libheader|Go-rcu}}
Timings: about 25 seconds for first 24 and 596 seconds for first 29.
<syntaxhighlight lang="go">package main
import (
"fmt"
big "github.com/ncw/gmp"
"rcu"
)
func main() {
const limit = 29
count := 0
p := 1
one := big.NewInt(1)
three := big.NewInt(3)
w := new(big.Int)
for count < limit {
for {
p += 2
if rcu.IsPrime(p) {
break
}
}
w.SetUint64(1)
w.Lsh(w, uint(p))
w.Add(w, one)
w.Quo(w, three)
if w.ProbablyPrime(15) {
count++
ws := w.String()
le := len(ws)
if le >= 34 {
ws = ws[0:15] + "..." + ws[le-15:]
}
fmt.Printf("%5d: %s", p, ws)
if le >= 34 {
fmt.Printf(" (%d digits)", le)
}
println()
}
}
}</syntaxhighlight>
{{out}}
<pre>
Same as C entry
</pre>
=={{header|J}}==
<syntaxhighlight lang=J>
11 683
13 2731
17 43691
19 174763
23 2796203
31 715827883
43 2932031007403
61 768614336404564651
79 201487636602438195784363
101 845100400152152934331135470251
127 56713727820156410577229101238628035243
167 62357403192785191176690552862561408838653121833643
191 1046183622564446793972631570534611069350392574077339085483
199 267823007376498379256993682056860433753700498963798805883563</syntaxhighlight>
Stretch goal (times are elapsed times, machine was a 1.2GHz intel i3):
<syntaxhighlight lang=J>{{
T0=. 6!:1''
f=. 3%~1+2x^]
c=. 0
echo 'count power seconds'
for_p. p:i.1e4 do.
if. 1 p: f p do.
c=. c+1
echo 5 6 8j3":c, p, (6!:1'')-T0
if. 24 <: c do. return. end.
end.
end.
}}_
count power seconds
1 3 0.004
2 5 0.006
3 7 0.008
4 11 0.013
5 13 0.018
6 17 0.021
7 19 0.025
8 23 0.028
9 31 0.030
10 43 0.033
11 61 0.035
12 79 0.039
13 101 0.044
14 127 0.052
15 167 0.062
16 191 0.075
17 199 0.089
18 313 0.120
19 347 0.167
20 701 0.436
21 1709 4.140
22 2617 16.035
23 3539 45.089
24 5807 181.280</syntaxhighlight>
=={{header|Julia}}==
* Requires Primes.jl
<syntaxhighlight lang=julia>
using Primes
function wagstaffpair(p::Integer)
isodd(p) || return (false, nothing)
isprime(p) || return (false, nothing)
m = (2^big(p) + 1) ÷ 3
isprime(m) || return (false, nothing)
return (true, m)
end
function findn_wagstaff_pairs(n_to_find::T) where T <: Integer
pairs = Tuple{T, BigInt}[]
count = 0
i = 2
while count < n_to_find
iswag, m = wagstaffpair(i)
iswag && push!(pairs, (i, m))
count += iswag
i += 1
end
return pairs
end
function println_wagstaff(pair; max_digit_display::Integer=20)
p, m = pair
mstr = string(m)
if length(mstr) > max_digit_display
umiddle = cld(max_digit_display, 2)
lmiddle = fld(max_digit_display, 2)
mstr = join((mstr[1:umiddle], "...", mstr[end-lmiddle+1:end],
" ($(length(mstr)) digits)"))
end
println("p = $p, m = $mstr")
end
foreach(println_wagstaff, findn_wagstaff_pairs(24))
</syntaxhighlight>
<pre>
p = 3, m = 3
p = 5, m = 11
p = 7, m = 43
p = 11, m = 683
p = 13, m = 2731
p = 17, m = 43691
p = 19, m = 174763
p = 23, m = 2796203
p = 31, m = 715827883
p = 43, m = 2932031007403
p = 61, m = 768614336404564651
p = 79, m = 2014876366...8195784363 (24 digits)
p = 101, m = 8451004001...1135470251 (30 digits)
p = 127, m = 5671372782...8628035243 (38 digits)
p = 167, m = 6235740319...3121833643 (50 digits)
p = 191, m = 1046183622...7339085483 (58 digits)
p = 199, m = 2678230073...8805883563 (60 digits)
p = 313, m = 5562466239...8130434731 (94 digits)
p = 347, m = 9556244233...1903606443 (104 digits)
p = 701, m = 3506757267...7823854251 (211 digits)
p = 1709, m = 9619252724...9070528171 (514 digits)
p = 2617, m = 2081504709...3435947691 (788 digits)
p = 3539, m = 7379609820...6486497963 (1065 digits)
p = 5807, m = 4018496237...6663568043 (1748 digits)
</pre>
=={{header|Java}}==
Line 125 ⟶ 815:
public class Main {
public static void main(String[] args) {
BigInteger d = new BigInteger("3"), a
int lmt = 25, sl, c = 0;
for (int i = 3; i < 5808; ) {
Line 166 ⟶ 856:
24 5807 40184962373..86663568043 1748 digits
</pre>
=={{header|Nim}}==
{{libheader|Nim-Integers}}
<syntaxhighlight lang=Nim>import std/strformat
import integers
func compress(str: string; size: int): string =
if str.len <= 2 * size: str
else: &"{str[0..<size]}...{str[^size..^1]} ({str.len} digits)"
echo "First 24 Wagstaff primes:"
let One = newInteger(1)
var count = 0
var p = 3
while count < 24:
if p.isPrime:
let n = (One shl p + 1) div 3
if n.isPrime:
inc count
echo &"{p:4}: {compress($n, 15)}"
inc p, 2
</syntaxhighlight>
{{out}}
<pre>First 24 Wagstaff primes:
3: 3
5: 11
7: 43
11: 683
13: 2731
17: 43691
19: 174763
23: 2796203
31: 715827883
43: 2932031007403
61: 768614336404564651
79: 201487636602438195784363
101: 845100400152152934331135470251
127: 567137278201564...101238628035243 (38 digits)
167: 623574031927851...838653121833643 (50 digits)
191: 104618362256444...574077339085483 (58 digits)
199: 267823007376498...963798805883563 (60 digits)
313: 556246623937737...099018130434731 (94 digits)
347: 955624423329196...712921903606443 (104 digits)
701: 350675726769891...862167823854251 (211 digits)
1709: 961925272487010...857299070528171 (514 digits)
2617: 208150470990263...847933435947691 (788 digits)
3539: 737960982013072...304686486497963 (1065 digits)
5807: 401849623734300...466686663568043 (1748 digits)
</pre>
=={{header|OCaml}}==
<syntaxhighlight lang="ocaml">let is_prime n =
let rec test x =
let q = n / x in x > q || x * q <> n && n mod (x + 2) <> 0 && test (x + 6)
in if n < 5 then n lor 1 = 3 else n land 1 <> 0 && n mod 3 <> 0 && test 5
let is_wagstaff n =
let w = succ (1 lsl n) / 3 in
if is_prime n && is_prime w then Some (n, w) else None
let () =
let show (p, w) = Printf.printf "%u -> %u%!\n" p w in
Seq.(ints 3 |> filter_map is_wagstaff |> take 11 |> iter show)</syntaxhighlight>
{{out}}
<pre>
3 -> 3
5 -> 11
7 -> 43
11 -> 683
13 -> 2731
17 -> 43691
19 -> 174763
23 -> 2796203
31 -> 715827883
43 -> 2932031007403
61 -> 768614336404564651
</pre>
=={{header|PARI/GP}}==
wprp.gp
<syntaxhighlight lang="c">\\ compiler: gp2c option: gp2c-run -g wprp.gp
/* wprp(p): input odd prime p > 5 . */
/* returns 1 if (2^p+1)/3 is a Wagstaff probable prime. */
wprp(p) = {
/* trial division up to a reasonable depth (time ratio tdiv/llr ≈ 0.2) */
my(l=log(p), ld=log(l));
forprimestep(q = 1, sqr(ld)^(l/log(2))\4, p+p,
if(Mod(2,q)^p == -1, return)
);
/* From R. & H. LIFCHITZ July 2000 */
/* see: http://primenumbers.net/Documents/TestNP.pdf */
/* if (2^p+1)/3 is prime ==> 25^2^(p-1) ≡ 25 (mod 2^p+1) */
/* Lucas-Lehmer-Riesel test with fast modular reduction. */
my(s=25, m=2^p-1);
for(i = 2, p,
s = sqr(s);
s = bitand(s,m) - s>>p
);
s==25
}; /* end wprp */
/* get exponents of Wagstaff prps in range [a,b] */
wprprun(a, b) = {
my(t=0, c=0, thr=default(nbthreads));
a = max(a,3);
gettime();
if(a <= 5,
if(a == 3, c++; p = 3; printf("#%d\tW%d\t%2dmin, %2d,%03d ms\n", c, p, t\60000%60, t\1000%60, t%1000));
c++; p = 5; printf("#%d\tW%d\t%2dmin, %2d,%03d ms\n", c, p, t\60000%60, t\1000%60, t%1000);
a = 7
);
parforprime(p= a, b, wprp(p), d, /* wprp(p) -> d copy from parallel world into real world. */
if(d,
t += gettime()\thr;
c++;
printf("#%d\tW%d\t%2dmin, %2d,%03d ms\n", c, p, t\60000%60, t\1000%60, t%1000)
)
)
}; /* end wprprun */
/* if running wprp as script */
\\ export(wprp);
wprprun(2, 42737)</syntaxhighlight>
{{out}}
Done on Intel(R) Core(TM) i5-8250U CPU @ 1.60GHz, 4 hyperthreaded cores.
<pre>#1 W3 0min, 0,000 ms
#2 W5 0min, 0,000 ms
#3 W7 0min, 0,000 ms
#4 W11 0min, 0,000 ms
#5 W13 0min, 0,000 ms
#6 W17 0min, 0,000 ms
#7 W19 0min, 0,000 ms
#8 W23 0min, 0,000 ms
#9 W31 0min, 0,000 ms
#10 W43 0min, 0,000 ms
#11 W61 0min, 0,000 ms
#12 W79 0min, 0,000 ms
#13 W101 0min, 0,000 ms
#14 W127 0min, 0,000 ms
#15 W167 0min, 0,000 ms
#16 W191 0min, 0,000 ms
#17 W199 0min, 0,000 ms
#18 W313 0min, 0,000 ms
#19 W347 0min, 0,000 ms
#20 W701 0min, 0,002 ms
#21 W1709 0min, 0,018 ms
#22 W2617 0min, 0,055 ms
#23 W3539 0min, 0,113 ms
#24 W5807 0min, 0,422 ms
#25 W10501 0min, 2,545 ms
#26 W10691 0min, 2,612 ms
#27 W11279 0min, 3,113 ms
#28 W12391 0min, 4,252 ms
#29 W14479 0min, 7,719 ms
#30 W42737 7min, 33,111 ms
? ##
*** last result: cpu time 1h, 8,965 ms, real time 7min, 33,449 ms.</pre>
=={{header|Perl}}==
Line 329 ⟶ 1,186:
24: 5807 => [1748 digit number]
</pre>
=={{header|Quackery}}==
<code>from</code>, <code>index</code>, <code>incr</code>, and <code>end</code> are defined at [[Loops/Increment loop index within loop body#Quackery]].
<code>prime</code> is defined at [[Miller–Rabin primality test#Quackery]].
<syntaxhighlight lang="Quackery"> [ bit 1+ 3 / ] is wagstaff ( n --> n )
[] 3 from
[ 2 incr
index prime while
index wagstaff prime while
index join
dup size 24 = if end ]
10 split swap
witheach
[ dup say "p = " echo
wagstaff
say ", w(p) = " echo cr ]
witheach
[ say "p = " echo cr ]
</syntaxhighlight>
{{out}}
<pre>p = 3, w(p) = 3
p = 5, w(p) = 11
p = 7, w(p) = 43
p = 11, w(p) = 683
p = 13, w(p) = 2731
p = 17, w(p) = 43691
p = 19, w(p) = 174763
p = 23, w(p) = 2796203
p = 31, w(p) = 715827883
p = 43, w(p) = 2932031007403
p = 61
p = 79
p = 101
p = 127
p = 167
p = 191
p = 199
p = 313
p = 347
p = 701
p = 1709
p = 2617
p = 3539
p = 5807</pre>
=={{header|Raku}}==
Line 382 ⟶ 1,289:
29: 14479 (179.87)
^C</pre>
=={{header|RPL}}==
{{works with|HP|49}}
≪ { } 2
'''WHILE''' OVER SIZE 10 < '''REPEAT'''
NEXTPRIME
2 OVER ^ 1 + 3 /
'''IF''' DUP ISPRIME? '''THEN'''
OVER →STR ":" + SWAP + ROT SWAP + SWAP
'''ELSE''' DROP '''END'''
'''END''' DROP
≫ '<span style="color:blue">WAGSTAFF</span>' STO
{{out}}
<pre>
1: { "3:3" "5:11" "7:43" "11:683" "13:2731" "17:43691" "19:174763" "23:2796203" "31:715827883" "43:2932031007403" }
</pre>
=={{header|Ruby}}==
For the first 12 the Ruby prime library is sufficiently fast. Here, the external GMP gem is used in order to do more.
<syntaxhighlight lang="ruby">require 'prime'
require 'gmp'
wagstaffs = Enumerator.new do |y|
odd_primes = Prime.each
odd_primes.next #skip 2
loop do
p = odd_primes.next
candidate = (2 ** p + 1)/3
y << [p, candidate] unless GMP::Z.new(candidate).probab_prime?.zero?
end
end
10.times{puts "%5d - %s" % wagstaffs.next}
14.times{puts "%5d" % wagstaffs.next.first}
</syntaxhighlight>
{{out}}
<pre> 3 - 3
5 - 11
7 - 43
11 - 683
13 - 2731
17 - 43691
19 - 174763
23 - 2796203
31 - 715827883
43 - 2932031007403
61
79
101
127
167
191
199
313
347
701
1709
2617
3539
5807
</pre>
=={{header|Wren}}==
Line 390 ⟶ 1,358:
Even so, as far as the stretch goal is concerned, takes around 25 seconds to find the next 14 Wagstaff primes but almost 10 minutes to find the next 19 on my machine (Core i7). I gave up after that.
<syntaxhighlight lang="
import "./gmp" for Mpz
import "./fmt" for Fmt
Line 472 ⟶ 1,440:
12391 (341 secs)
14479 (593 secs)
</pre>
=={{header|XPL0}}==
<syntaxhighlight lang "XPL0">func IsPrime(N); \Return 'true' if N is prime
real N; int I;
[if N <= 2. then return N = 2.;
if Mod(N, 2.) = 0. then \even\ return false;
for I:= 3 to fix(sqrt(N)) do
[if Mod(N, float(I)) = 0. then return false;
I:= I+1;
];
return true;
];
real P, Q; int C;
[P:= 2.; C:= 0;
Format(1, 0);
repeat if IsPrime(P) then
[Q:= Pow(2., P) + 1.;
if Mod(Q, 3.) = 0. and IsPrime(Q/3.) then
[Text(0, "(2^^");
RlOut(0, P);
Text(0, " + 1)/3 = ");
RlOut(0, Q/3.);
CrLf(0);
C:= C+1;
];
];
P:= P+1.;
until C >= 10;
]</syntaxhighlight>
{{out}}
<pre>
(2^3 + 1)/3 = 3
(2^5 + 1)/3 = 11
(2^7 + 1)/3 = 43
(2^11 + 1)/3 = 683
(2^13 + 1)/3 = 2731
(2^17 + 1)/3 = 43691
(2^19 + 1)/3 = 174763
(2^23 + 1)/3 = 2796203
(2^31 + 1)/3 = 715827883
(2^43 + 1)/3 = 2932031007403
</pre>
| |||