Rosetta Code

Pierpont primes: Difference between revisions

Page Discussion

Pierpont primes (view source)

Revision as of 11:47, 25 January 2024

78,047 bytes added ,  3 months ago
m
→‎{{header|Wren}}: Minor tidy
m (Fix Perl 6 -> Raku links)
m (→‎{{header|Wren}}: Minor tidy)
 
(48 intermediate revisions by 26 users not shown)
Line 28:
 
<br>
 
=={{header|ALGOL 68}}==
As in the second Go sample, generates the 3-smooth number sequence to find Pierpoint Prime candidates. Uses a sliding window on the sequence to avoid having to keep it all in memory.
{{libheader|ALGOL 68-primes}}
<syntaxhighlight lang="algol68">BEGIN # find some pierpoint primes of the first kind (2^u3^v + 1) #
# and second kind (2^u3^v - 1) #
# construct a sieve of primes up to 10 000 #
PR read "primes.incl.a68" PR
[]BOOL prime = PRIMESIEVE 10 000;
# returns the minimum of a and b #
PROC min = ( INT a, b )INT: IF a < b THEN a ELSE b FI;
# returns TRUE if p is prime, FALSE otherwise #
PROC is prime = ( INT p )BOOL:
IF NOT ODD p
THEN p = 2
ELIF p <= UPB prime
THEN prime[ p ] # small enough to use the sieve #
ELSE # use trial division by the sieved primes #
BOOL probably prime := TRUE;
FOR d FROM 3 BY 2 WHILE d * d <= p AND probably prime DO
IF prime[ d ] THEN probably prime := p MOD d /= 0 FI
OD;
probably prime
FI; # is prime #
# sets the elements of pp1 to the first UPB pp1 Pierpoint primes of the #
# first kind and pp2 to the first UPB pp2 Pierpoint primes of the #
# second kind #
PROC find pierpoint = ( REF[]INT pp1, pp2 )VOID:
BEGIN
# saved 3-smooth values #
# - only the currently active part of the seuence is kept #
INT three smooth limit = 33;
[ 1 : three smooth limit * 3 ]INT s3s; FOR i TO UPB s3s DO s3s[ i ] := 0 OD;
INT pp1 pos := LWB pp1 - 1, pp2 pos := LWB pp2 - 1;
INT pp1 max = UPB pp1;
INT pp2 max = UPB pp2;
INT p2, p3, last2, last3;
INT s pos := s3s[ 1 ] := last2 := last3 := p2 := p3 := 1;
FOR n FROM 2 WHILE pp1 pos < pp1 max OR pp2 pos < pp2 max DO
# the next 3-smooth number is the lowest of the next #
# multiples of 2 and 3 #
INT m = min( p2, p3 );
IF pp1 pos < pp1 max THEN
IF INT first = m + 1;
is prime( first )
THEN # have a Pierpoint prime of the first kind #
pp1[ pp1 pos +:= 1 ] := first
FI
FI;
IF pp2 pos < pp2 max THEN
IF INT second = m - 1;
is prime( second )
THEN # have a Pierpoint prime of the second kind #
pp2[ pp2 pos +:= 1 ] := second
FI
FI;
s3s[ s pos +:= 1 ] := m;
IF m = p2 THEN p2 := 2 * s3s[ last2 +:= 1 ] FI;
IF m = p3 THEN p3 := 3 * s3s[ last3 +:= 1 ] FI;
INT last min = IF last2 > last3 THEN last3 ELSE last2 FI;
IF last min > three smooth limit THEN
# shuffle the sequence down, over the now unused bits #
INT new pos := 0;
FOR pos FROM last min TO s pos DO
s3s[ new pos +:= 1 ] := s3s[ pos ]
OD;
INT diff := last min - 1;
last2 -:= diff;
last3 -:= diff;
s pos -:= diff
FI
OD
END; # find pierpoint #
# prints a table of Pierpoint primes of the specified kind #
PROC print pierpoint = ( []INT primes, STRING kind )VOID:
BEGIN
print( ( "The first "
, whole( ( UPB primes + 1 ) - LWB primes, 0 )
, " Pierpoint primes of the "
, kind
, " kind:"
, newline
)
);
INT p count := 0;
FOR i FROM LWB primes TO UPB primes DO
print( ( " ", whole( primes[ i ], -8 ) ) );
IF ( p count +:= 1 ) MOD 10 = 0 THEN print( ( newline ) ) FI
OD;
print( ( newline ) )
END; # print pierpoint #
# find the first 50 Pierpoint primes of the first and second kinds #
INT max pierpoint = 50;
[ 1 : max pierpoint ]INT pfirst;
[ 1 : max pierpoint ]INT psecond;
find pierpoint( pfirst, psecond );
# show the Pierpoint primes #
print pierpoint( pfirst, "First" );
print pierpoint( psecond, "Second" )
END</syntaxhighlight>
{{out}}
<pre>
The first 50 Pierpoint primes of the First kind:
2 3 5 7 13 17 19 37 73 97
109 163 193 257 433 487 577 769 1153 1297
1459 2593 2917 3457 3889 10369 12289 17497 18433 39367
52489 65537 139969 147457 209953 331777 472393 629857 746497 786433
839809 995329 1179649 1492993 1769473 1990657 2654209 5038849 5308417 8503057
 
The first 50 Pierpoint primes of the Second kind:
2 3 5 7 11 17 23 31 47 53
71 107 127 191 383 431 647 863 971 1151
2591 4373 6143 6911 8191 8747 13121 15551 23327 27647
62207 73727 131071 139967 165887 294911 314927 442367 472391 497663
524287 786431 995327 1062881 2519423 10616831 17915903 18874367 25509167 30233087
</pre>
 
=={{header|C}}==
{{trans|D}}
<syntaxhighlight lang="c">#include <stdbool.h>
#include <stdint.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
 
const int PRIMES[] = {
2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47,
53, 59, 61, 67, 71, 73, 79, 83, 89, 97, 101, 103, 107, 109, 113,
127, 131, 137, 139, 149, 151, 157, 163, 167, 173, 179, 181, 191, 193, 197,
199, 211, 223, 227, 229, 233, 239, 241, 251, 257, 263, 269, 271, 277, 281,
283, 293, 307, 311, 313, 317, 331, 337, 347, 349, 353, 359, 367, 373, 379,
383, 389, 397, 401, 409, 419, 421, 431, 433, 439, 443, 449, 457, 461, 463,
467, 479, 487, 491, 499, 503, 509, 521, 523, 541, 547, 557, 563, 569, 571,
577, 587, 593, 599, 601, 607, 613, 617, 619, 631, 641, 643, 647, 653, 659,
661, 673, 677, 683, 691, 701, 709, 719, 727, 733, 739, 743, 751, 757, 761,
769, 773, 787, 797, 809, 811, 821, 823, 827, 829, 839, 853, 857, 859, 863,
877, 881, 883, 887, 907, 911, 919, 929, 937, 941, 947, 953, 967, 971, 977,
};
#define PRIME_SIZE (sizeof(PRIMES) / sizeof(int))
 
bool isPrime(const int n) {
int i;
 
if (n < 2) {
return false;
}
 
for (i = 0; i < PRIME_SIZE; i++) {
if (n == PRIMES[i]) {
return true;
}
if (n % PRIMES[i] == 0) {
return false;
}
}
 
if (n < PRIMES[PRIME_SIZE - 1] * PRIMES[PRIME_SIZE - 1]) {
return true;
}
 
i = PRIMES[PRIME_SIZE - 1]+2;
while (i * i < n) {
if (n % i == 0) {
return false;
}
i += 2;
}
 
return true;
}
 
#define N 50
int p[2][50];
void pierpont() {
int64_t s[8 * N];
int count = 0;
int count1 = 1;
int count2 = 0;
int i2 = 0;
int i3 = 0;
int k = 1;
int64_t n2, n3, t;
int64_t *sp = &s[1];
 
memset(p[0], 0, N * sizeof(int));
memset(p[1], 0, N * sizeof(int));
p[0][0] = 2;
s[0] = 1;
 
while (count < N) {
n2 = s[i2] * 2;
n3 = s[i3] * 3;
if (n2 < n3) {
t = n2;
i2++;
} else {
t = n3;
i3++;
}
if (t > s[k - 1]) {
*sp++ = t;
k++;
 
t++;
if (count1 < N && isPrime(t)) {
p[0][count1] = t;
count1++;
}
 
t -= 2;
if (count2 < N && isPrime(t)) {
p[1][count2] = t;
count2++;
}
 
count = min(count1, count2);
}
}
}
 
int main() {
int i;
 
pierpont();
 
printf("First 50 Pierpont primes of the first kind:\n");
for (i = 0; i < N; i++) {
printf("%8d ", p[0][i]);
if ((i - 9) % 10 == 0) {
printf("\n");
}
}
printf("\n");
 
printf("First 50 Pierpont primes of the second kind:\n");
for (i = 0; i < N; i++) {
printf("%8d ", p[1][i]);
if ((i - 9) % 10 == 0) {
printf("\n");
}
}
printf("\n");
}</syntaxhighlight>
{{out}}
<pre>First 50 Pierpont primes of the first kind:
2 3 5 7 13 17 19 37 73 97
109 163 193 257 433 487 577 769 1153 1297
1459 2593 2917 3457 3889 10369 12289 17497 18433 39367
52489 65537 139969 147457 209953 331777 472393 629857 746497 786433
839809 995329 1179649 1492993 1769473 1990657 2654209 5038849 5308417 8503057
 
First 50 Pierpont primes of the second kind:
2 3 5 7 11 17 23 31 47 53
71 107 127 191 383 431 647 863 971 1151
2591 4373 6143 6911 8191 8747 13121 15551 23327 27647
62207 73727 131071 139967 165887 294911 314927 442367 472391 497663
524287 786431 995327 1062881 2519423 10616831 17915903 18874367 25509167 30233087</pre>
 
=={{header|C sharp|C#}}==
<langsyntaxhighlight lang="csharp">using System;
using System.Collections.Generic;
using System.Numerics;
Line 202 ⟶ 459:
}
}
}</langsyntaxhighlight>
{{out}}
<pre>First 50 Pierpont primes of the first kind:
Line 221 ⟶ 478:
250th Pierpont prime of the second kind: 4111131172000956525894875083702271</pre>
 
=={{header|C++}}==
{{libheader|GMP}}
<syntaxhighlight lang="cpp">#include <cassert>
#include <algorithm>
#include <iomanip>
#include <iostream>
#include <vector>
#include <gmpxx.h>
 
using big_int = mpz_class;
 
bool is_prime(const big_int& n) {
return mpz_probab_prime_p(n.get_mpz_t(), 25);
}
 
template <typename integer>
class n_smooth_generator {
public:
explicit n_smooth_generator(size_t n);
integer next();
size_t size() const {
return results_.size();
}
private:
std::vector<size_t> primes_;
std::vector<size_t> index_;
std::vector<integer> results_;
std::vector<integer> queue_;
};
 
template <typename integer>
n_smooth_generator<integer>::n_smooth_generator(size_t n) {
for (size_t p = 2; p <= n; ++p) {
if (is_prime(p)) {
primes_.push_back(p);
queue_.push_back(p);
}
}
index_.assign(primes_.size(), 0);
results_.push_back(1);
}
 
template <typename integer>
integer n_smooth_generator<integer>::next() {
integer last = results_.back();
for (size_t i = 0, n = primes_.size(); i < n; ++i) {
if (queue_[i] == last)
queue_[i] = results_[++index_[i]] * primes_[i];
}
results_.push_back(*min_element(queue_.begin(), queue_.end()));
return last;
}
 
void print_vector(const std::vector<big_int>& numbers) {
for (size_t i = 0, n = numbers.size(); i < n; ++i) {
std::cout << std::setw(9) << numbers[i];
if ((i + 1) % 10 == 0)
std::cout << '\n';
}
std::cout << '\n';
}
 
int main() {
const int max = 250;
std::vector<big_int> first, second;
int count1 = 0;
int count2 = 0;
n_smooth_generator<big_int> smooth_gen(3);
big_int p1, p2;
 
while (count1 < max || count2 < max) {
big_int n = smooth_gen.next();
if (count1 < max && is_prime(n + 1)) {
p1 = n + 1;
if (first.size() < 50)
first.push_back(p1);
++count1;
}
if (count2 < max && is_prime(n - 1)) {
p2 = n - 1;
if (second.size() < 50)
second.push_back(p2);
++count2;
}
}
std::cout << "First 50 Pierpont primes of the first kind:\n";
print_vector(first);
std::cout << "First 50 Pierpont primes of the second kind:\n";
print_vector(second);
std::cout << "250th Pierpont prime of the first kind: " << p1 << '\n';
std::cout << "250th Pierpont prime of the second kind: " << p2 << '\n';
return 0;
}</syntaxhighlight>
 
{{out}}
<pre>
First 50 Pierpont primes of the first kind:
2 3 5 7 13 17 19 37 73 97
109 163 193 257 433 487 577 769 1153 1297
1459 2593 2917 3457 3889 10369 12289 17497 18433 39367
52489 65537 139969 147457 209953 331777 472393 629857 746497 786433
839809 995329 1179649 1492993 1769473 1990657 2654209 5038849 5308417 8503057
 
First 50 Pierpont primes of the second kind:
2 3 5 7 11 17 23 31 47 53
71 107 127 191 383 431 647 863 971 1151
2591 4373 6143 6911 8191 8747 13121 15551 23327 27647
62207 73727 131071 139967 165887 294911 314927 442367 472391 497663
524287 786431 995327 1062881 2519423 10616831 17915903 18874367 25509167 30233087
 
250th Pierpont prime of the first kind: 62518864539857068333550694039553
250th Pierpont prime of the second kind: 4111131172000956525894875083702271
</pre>
 
=={{header|D}}==
{{trans|C#}}
<syntaxhighlight lang="d">import std.algorithm;
import std.bigint;
import std.random;
import std.stdio;
 
immutable PRIMES = [
2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47,
53, 59, 61, 67, 71, 73, 79, 83, 89, 97, 101, 103, 107, 109, 113,
127, 131, 137, 139, 149, 151, 157, 163, 167, 173, 179, 181, 191, 193, 197,
199, 211, 223, 227, 229, 233, 239, 241, 251, 257, 263, 269, 271, 277, 281,
283, 293, 307, 311, 313, 317, 331, 337, 347, 349, 353, 359, 367, 373, 379,
383, 389, 397, 401, 409, 419, 421, 431, 433, 439, 443, 449, 457, 461, 463,
467, 479, 487, 491, 499, 503, 509, 521, 523, 541, 547, 557, 563, 569, 571,
577, 587, 593, 599, 601, 607, 613, 617, 619, 631, 641, 643, 647, 653, 659,
661, 673, 677, 683, 691, 701, 709, 719, 727, 733, 739, 743, 751, 757, 761,
769, 773, 787, 797, 809, 811, 821, 823, 827, 829, 839, 853, 857, 859, 863,
877, 881, 883, 887, 907, 911, 919, 929, 937, 941, 947, 953, 967, 971, 977,
];
 
BigInt getRandom(BigInt min, BigInt max) {
auto r = max - min;
auto hs = r.toHex;
 
BigInt result;
do {
string t = "0x";
for (int i = 0; i < hs.length; i++) {
t ~= "0123456789abcdef"[uniform(0, 16)];
}
result = BigInt(t) + min;
} while (result < min || max <= result);
return result;
}
 
//Modified from https://rosettacode.org/wiki/Miller-Rabin_primality_test#Python
bool isProbablePrime(BigInt n) {
if (n == 0 || n == 1) {
return false;
}
 
bool check(BigInt num) {
foreach (prime; PRIMES) {
if (num == prime) {
return true;
}
if (num % prime == 0) {
return false;
}
if (prime * prime > num) {
return true;
}
}
return true;
}
 
if (check(n)) {
auto large = PRIMES[$ - 1];
if (n <= large) {
return true;
}
}
 
int s = 0;
auto d = n - 1;
while ((d & 1) == 0) {
d >>= 1;
s++;
}
 
bool trialComposite(BigInt a) {
if (powmod(a, d, n) == 1) {
return false;
}
for (int i = 0; i < s; i++) {
auto t = BigInt(2) ^^ i;
if (powmod(a, t * d, n) == n - 1) {
return false;
}
}
return true;
}
 
for (int i = 0; i < 8; i++) {
auto a = getRandom(BigInt(2), n);
if (trialComposite(a)) {
return false;
}
}
return true;
}
 
BigInt[][] pierpont(int n) {
BigInt[][] p = [[], []];
for (int i = 0; i < n; i++) {
p[0] ~= BigInt(0);
p[1] ~= BigInt(0);
}
p[0][0] = 2;
 
int count = 0;
int count1 = 1;
int count2 = 0;
BigInt[] s = [BigInt(1)];
int i2 = 0;
int i3 = 0;
int k = 1;
BigInt n2, n3, t;
 
while (count < n) {
n2 = s[i2] * 2;
n3 = s[i3] * 3;
if (n2 < n3) {
t = n2;
i2++;
} else {
t = n3;
i3++;
}
if (t > s[k - 1]) {
s ~= t;
k++;
 
t++;
if (count1 < n && t.isProbablePrime()) {
p[0][count1] = t;
count1++;
}
 
t -= 2;
if (count2 < n && t.isProbablePrime()) {
p[1][count2] = t;
count2++;
}
 
count = min(count1, count2);
}
}
 
return p;
}
 
void main() {
auto p = pierpont(250);
 
writeln("First 50 Pierpont primes of the first kind:");
for (int i = 0; i < 50; i++) {
writef("%8d ", p[0][i]);
if ((i - 9) % 10 == 0) {
writeln;
}
}
writeln;
 
writeln("First 50 Pierpont primes of the first kind:");
for (int i = 0; i < 50; i++) {
writef("%8d ", p[1][i]);
if ((i - 9) % 10 == 0) {
writeln;
}
}
writeln;
 
writefln("%dth Pierpont prime of the first kind: %d", p[0].length, p[0][$ - 1]);
writefln("%dth Pierpont prime of the second kind: %d", p[1].length, p[1][$ - 1]);
}</syntaxhighlight>
{{out}}
<pre>First 50 Pierpont primes of the first kind:
2 3 5 7 13 17 19 37 73 97
109 163 193 257 433 487 577 769 1153 1297
1459 2593 2917 3457 3889 10369 12289 17497 18433 39367
52489 65537 139969 147457 209953 331777 472393 629857 746497 786433
839809 995329 1179649 1492993 1769473 1990657 2654209 5038849 5308417 8503057
 
First 50 Pierpont primes of the first kind:
2 3 5 7 11 17 23 31 47 53
71 107 127 191 383 431 647 863 971 1151
2591 4373 6143 6911 8191 8747 13121 15551 23327 27647
62207 73727 131071 139967 165887 294911 314927 442367 472391 497663
524287 786431 995327 1062881 2519423 10616831 17915903 18874367 25509167 30233087
 
250th Pierpont prime of the first kind: 62518864539857068333550694039553
250th Pierpont prime of the second kind: 4111131172000956525894875083702271</pre>
=={{header|Delphi}}==
{{libheader| System.SysUtils}}
{{libheader| System.Math}}
{{libheader| System.StrUtils}}
{{libheader| System.Generics.Collections}}
{{libheader| System.Generics.Defaults}}
{{libheader| Velthuis.BigIntegers}}
{{libheader| Velthuis.BigIntegers.Primes}}
{{Trans|Go}}
Thanks Rudy Velthuis for the [https://github.com/rvelthuis/DelphiBigNumbers Velthuis.BigIntegers.Primes and Velthuis.BigIntegers] library.<br>
<syntaxhighlight lang="delphi">
program Pierpont_primes;
 
{$APPTYPE CONSOLE}
 
uses
System.SysUtils,
System.Math,
System.StrUtils,
System.Generics.Collections,
System.Generics.Defaults,
Velthuis.BigIntegers,
Velthuis.BigIntegers.Primes;
 
function Pierpont(ulim, vlim: Integer; first: boolean): TArray<BigInteger>;
begin
var p: BigInteger := 0;
var p2: BigInteger := 1;
var p3: BigInteger := 1;
 
for var v := 0 to vlim - 1 do
begin
for var u := 0 to ulim - 1 do
begin
p := p2 * p3;
if first then
p := p + 1
else
p := p - 1;
if IsProbablePrime(p, 10) then
begin
SetLength(result, Length(result) + 1);
result[High(result)] := BigInteger(p);
end;
p2 := p2 * 2;
end;
p3 := p3 * 3;
p2 := 1;
end;
 
TArray.sort<BigInteger>(Result, TComparer<BigInteger>.Construct(
function(const Left, Right: BigInteger): Integer
begin
Result := BigInteger.Compare(Left, Right);
end));
end;
 
begin
 
writeln('First 50 Pierpont primes of the first kind:');
var pp := Pierpont(120, 80, True);
for var i := 0 to 49 do
begin
write(pp[i].ToString: 8, ' ');
if ((i - 9) mod 10) = 0 then
writeln;
end;
 
writeln('First 50 Pierpont primes of the second kind:');
var pp2 := Pierpont(120, 80, False);
for var i := 0 to 49 do
begin
write(pp2[i].ToString: 8, ' ');
if ((i - 9) mod 10) = 0 then
writeln;
end;
 
Writeln('250th Pierpont prime of the first kind:', pp[249].ToString);
Writeln('250th Pierpont prime of the second kind:', pp2[249].ToString);
 
readln;
end.</syntaxhighlight>
=={{header|F_Sharp|F#}}==
This task uses [[Extensible_prime_generator#The_functions|Extensible Prime Generator (F#)]].<br>
<syntaxhighlight lang="fsharp">
// Pierpont primes . Nigel Galloway: March 19th., 2021
let fN g=let mutable g=g in ((fun()->g),fun()->g<-g+g;())
let fG y=let rec fG n g=seq{match g|>List.minBy(fun(n,_)->n()) with (f,s) when f()=n->yield f()+y; s(); yield! fG(n*3)(fN(n*3)::g)
|(f,s) ->yield f()+y; s(); yield! fG n g}
seq{yield! fG 1 [fN 1]}|>Seq.filter isPrime
let pierpontT1,pierpontT2=fG 1,fG -1
 
pierpontT1|>Seq.take 50|>Seq.iter(printf "%d "); printfn ""
pierpontT2|>Seq.take 50|>Seq.iter(printf "%d "); printfn ""
</syntaxhighlight>
{{out}}
<pre>
2 3 5 7 13 17 19 37 73 97 109 163 193 257 433 487 577 769 1153 1297 1459 2593 2917 3457 3889 10369 12289 17497 18433 39367 52489 65537 139969 147457 209953 331777 472393 629857 746497 786433 839809 995329 1179649 1492993 1769473 1990657 2654209 5038849 5308417 8503057
2 3 5 7 11 17 23 31 47 53 71 107 127 191 383 431 647 863 971 1151 2591 4373 6143 6911 8191 8747 13121 15551 23327 27647 62207 73727 131071 139967 165887 294911 314927 442367 472391 497663 524287 786431 995327 1062881 2519423 10616831 17915903 18874367 25509167 30233087
</pre>
=={{header|Factor}}==
<langsyntaxhighlight lang="factor">USING: fry grouping io kernel locals make math math.functions
math.primes prettyprint sequences sorting ;
 
Line 250 ⟶ 907:
"250th Pierpont prime of the second kind: " write
249 second nth .
]</langsyntaxhighlight>
{{out}}
<pre>
Line 270 ⟶ 927:
 
250th Pierpont prime of the second kind: 4111131172000956525894875083702271
</pre>
 
=={{header|FreeBASIC}}==
<syntaxhighlight lang="freebasic">#define NPP 50
 
Function isPrime(Byval n As Ulongint) As Boolean
If n < 2 Then Return false
If n = 2 Then Return true
If n Mod 2 = 0 Then Return false
For i As Uinteger = 3 To Int(Sqr(n))+1 Step 2
If n Mod i = 0 Then Return false
Next i
Return true
End Function
 
Function is_23(Byval n As Uinteger) As Boolean
While n Mod 2 = 0
n /= 2
Wend
While n Mod 3 = 0
n /= 3
Wend
Return Iif(n=1, true, false)
End Function
 
Function isPierpont(n As Uinteger) As Uinteger
If Not isPrime(n) Then Return 0 'not prime
Dim As Uinteger p1 = is_23(n+1), p2 = is_23(n-1)
If p1 And p2 Then Return 3 'pierpont prime of both kinds
If p1 Then Return 1 'pierpont prime of the 1st kind
If p2 Then Return 2 'pierpont prime of the 2nd kind
Return 0 'prime, but not pierpont
End Function
 
Dim As Uinteger pier(1 To 2, 1 To NPP), np(1 To 2) = {0, 0}
Dim As Uinteger x = 1, j
While np(1) <= NPP Or np(2) <= NPP
x += 1
j = isPierpont(x)
If j > 0 Then
If j Mod 2 = 1 Then
np(1) += 1
If np(1) <= NPP Then pier(1, np(1)) = x
End If
If j > 1 Then
np(2) += 1
If np(2) <= NPP Then pier(2, np(2)) = x
End If
End If
Wend
 
Print "First 50 Pierpoint primes of the first kind:"
For j = 1 To NPP
Print Using " ########"; pier(2, j);
If j Mod 10 = 0 Then Print
Next j
Print !"\nFirst 50 Pierpoint primes of the secod kind:"
For j = 1 To NPP
Print Using " ########"; pier(1, j);
If j Mod 10 = 0 Then Print
Next j
</syntaxhighlight>
{{out}}
<pre>
First 50 Pierpont primes of the first kind:
2 3 5 7 13 17 19 37 73 97
109 163 193 257 433 487 577 769 1153 1297
1459 2593 2917 3457 3889 10369 12289 17497 18433 39367
52489 65537 139969 147457 209953 331777 472393 629857 746497 786433
839809 995329 1179649 1492993 1769473 1990657 2654209 5038849 5308417 8503057
 
First 50 Pierpont primes of the second kind:
2 3 5 7 11 17 23 31 47 53
71 107 127 191 383 431 647 863 971 1151
2591 4373 6143 6911 8191 8747 13121 15551 23327 27647
62207 73727 131071 139967 165887 294911 314927 442367 472391 497663
524287 786431 995327 1062881 2519423 10616831 17915903 18874367 25509167 30233087
</pre>
 
Line 281 ⟶ 1,015:
 
However, in order to be sure that the first 250 Pierpont primes will be generated, it is necessary to choose the loop sizes to produce somewhat more than this and then sort them into order.
<langsyntaxhighlight lang="go">package main
 
import (
Line 342 ⟶ 1,076:
fmt.Println("\n250th Pierpont prime of the first kind:", pp[249])
fmt.Println("\n250th Pierpont prime of the second kind:", pp2[249])
}</langsyntaxhighlight>
 
{{out}}
Line 374 ⟶ 1,108:
These timings are for my Celeron @1.6GHz and should therefore be much faster on a more modern machine.
 
<langsyntaxhighlight lang="go">package main
 
import (
Line 466 ⟶ 1,200:
elapsed := time.Now().Sub(start)
fmt.Printf("\nTook %s\n", elapsed)
}</langsyntaxhighlight>
 
{{out}}
Line 498 ⟶ 1,232:
Took 1m40.781726122s
</pre>
 
=={{header|Haskell}}==
Uses arithmoi Library: https://hackage.haskell.org/package/arithmoi-0.11.0.0 for prime generation and prime testing.<br>
n-smooth generation function based on [[Hamming_numbers#Avoiding_generation_of_duplicates]]
<syntaxhighlight lang="haskell">import Control.Monad (guard)
import Data.List (intercalate)
import Data.List.Split (chunksOf)
import Math.NumberTheory.Primes (Prime, unPrime, nextPrime)
import Math.NumberTheory.Primes.Testing (isPrime)
import Text.Printf (printf)
 
data PierPointKind = First | Second
 
merge :: Ord a => [a] -> [a] -> [a]
merge [] b = b
merge a@(x:xs) b@(y:ys) | x < y = x : merge xs b
| otherwise = y : merge a ys
 
nSmooth :: Integer -> [Integer]
nSmooth p = 1 : foldr u [] factors
where
factors = takeWhile (<=p) primes
primes = map unPrime [nextPrime 1..]
u n s = r
where
r = merge s (map (n*) (1:r))
 
pierpoints :: PierPointKind -> [Integer]
pierpoints k = do
n <- nSmooth 3
let x = case k of First -> succ n
Second -> pred n
guard (isPrime x) >> [x]
 
main :: IO ()
main = do
printf "\nFirst 50 Pierpont primes of the first kind:\n"
mapM_ (\row -> mapM_ (printf "%12s" . commas) row >> printf "\n") (rows $ pierpoints First)
printf "\nFirst 50 Pierpont primes of the second kind:\n"
mapM_ (\row -> mapM_ (printf "%12s" . commas) row >> printf "\n") (rows $ pierpoints Second)
printf "\n250th Pierpont prime of the first kind: %s\n" (commas $ pierpoints First !! 249)
printf "\n250th Pierpont prime of the second kind: %s\n\n" (commas $ pierpoints Second !! 249)
where
rows = chunksOf 10 . take 50
commas = reverse . intercalate "," . chunksOf 3 . reverse . show</syntaxhighlight>
{{out}}
<pre>
First 50 Pierpont primes of the first kind:
2 3 5 7 13 17 19 37 73 97
109 163 193 257 433 487 577 769 1,153 1,297
1,459 2,593 2,917 3,457 3,889 10,369 12,289 17,497 18,433 39,367
52,489 65,537 139,969 147,457 209,953 331,777 472,393 629,857 746,497 786,433
839,809 995,329 1,179,649 1,492,993 1,769,473 1,990,657 2,654,209 5,038,849 5,308,417 8,503,057
 
First 50 Pierpont primes of the second kind:
2 3 5 7 11 17 23 31 47 53
71 107 127 191 383 431 647 863 971 1,151
2,591 4,373 6,143 6,911 8,191 8,747 13,121 15,551 23,327 27,647
62,207 73,727 131,071 139,967 165,887 294,911 314,927 442,367 472,391 497,663
524,287 786,431 995,327 1,062,881 2,519,423 10,616,831 17,915,903 18,874,367 25,509,167 30,233,087
 
250th Pierpont prime of the first kind: 62,518,864,539,857,068,333,550,694,039,553
 
250th Pierpont prime of the second kind: 4,111,131,172,000,956,525,894,875,083,702,271
 
./pierpoints 0.04s user 0.01s system 20% cpu 0.215 total
</pre>
 
=={{header|J}}==
 
Implementation (first kind first):
 
<syntaxhighlight lang="j"> 5 10$(#~ 1 p:])1+/:~,*//2 3x^/i.20
2 3 5 7 13 17 19 37 73 97
109 163 193 257 433 487 577 769 1153 1297
1459 2593 2917 3457 3889 10369 12289 17497 18433 39367
52489 65537 139969 147457 209953 331777 472393 629857 746497 786433
839809 995329 1179649 1492993 1769473 1990657 2654209 5038849 5308417 8503057
5 10$(#~ 1 p:])_1+/:~,*//2 3x^/i.20
2 3 5 7 11 17 23 31 47 53
71 107 127 191 383 431 647 863 971 1151
2591 4373 6143 6911 8191 8747 13121 15551 23327 27647
62207 73727 131071 139967 165887 294911 314927 442367 472391 497663
524287 786431 995327 1062881 2519423 10616831 17915903 25509167 30233087 57395627</syntaxhighlight>
 
and
 
<syntaxhighlight lang="j"> 249{ (#~ 1 p:])1+/:~,*//2 3x^/i.112
62518864539857068333550694039553
249{ (#~ 1 p:])_1+/:~,*//2 3x^/i.112
4111131172000956525894875083702271</syntaxhighlight>
 
=={{header|Java}}==
<langsyntaxhighlight lang="java">
import java.math.BigInteger;
import java.text.NumberFormat;
Line 566 ⟶ 1,391:
}
</syntaxhighlight>
</lang>
 
{{out}}
Line 587 ⟶ 1,412:
 
250th Pierpont prime of the second kind: 4,111,131,172,000,956,525,894,875,083,702,271
</pre>
 
=={{header|jq}}==
{{works with|jq}}
'''Works with gojq, the Go implementation of jq'''
 
Since we do not know how often 2^u * 3^v ± 1 is prime, we use an
adaptive approach based on adjusting the upper bounds for u and v
that we need consider. The idea is to avoid any hand-waving about the
correctness of the solution. Some "debug" statements that are
informative about the adaptive steps have been retained as comments.
 
The standard (C-based) implementation of jq does not support very large
integers, and the Go implementation, gojq, requires too much memory
to solve the stretch problem, so the code for the stretch task is
shown, but not executed.
 
'''Preliminaries'''
<syntaxhighlight lang="jq">def is_prime:
. as $n
| if ($n < 2) then false
elif ($n % 2 == 0) then $n == 2
elif ($n % 3 == 0) then $n == 3
elif ($n % 5 == 0) then $n == 5
elif ($n % 7 == 0) then $n == 7
elif ($n % 11 == 0) then $n == 11
elif ($n % 13 == 0) then $n == 13
elif ($n % 17 == 0) then $n == 17
elif ($n % 19 == 0) then $n == 19
else {i:23}
| until( (.i * .i) > $n or ($n % .i == 0); .i += 2)
| .i * .i > $n
end;
 
# pretty-printing
def lpad($len): tostring | ($len - length) as $l | (" " * $l)[:$l] + .;
 
def nwise($n):
def n: if length <= $n then . else .[0:$n] , (.[$n:] | n) end;
n;
 
def table($ncols; $colwidth):
nwise($ncols) | map(lpad($colwidth)) | join(" ");</syntaxhighlight>
'''Pierpont primes'''
<syntaxhighlight lang="jq"># Input: the target number of primes of the form p(u,v) == 2^u * 3^v ± 1.
# The main idea is to let u run from 0:m and v run from 0:n where n ~~ 0.63 * m.
# Initially we start with plausible values for m and n ($maxm and $maxn respectively),
# and then check whether these have been chosen conservatively enough.
#
def pierponts($firstkind):
. as $N
| (if $firstkind then 1 else -1 end) as $incdec
 
# Input: [$maxm, $maxn]
# Output: an array of objects of the form {p, m, n}
# where .p is prime and .m and .n are the corresponding powers of 2 and 3
| def pp:
. as [$maxm, $maxn]
| [ ({p2:1, m:0},
(foreach range(0; $maxm) as $m (1; . * 2; {p2: ., m: ($m + 1)}))) as $a
| ({p3:1, n:0},
(foreach range(0; $maxn) as $n (1; . * 3; {p3: ., n: ($n + 1)}))) as $b
| {p: ($a.p2 * $b.p3 + $incdec), m: $a.m, n: $b.n}
| select(.p|is_prime) ]
| unique_by(.p)
# | (length|debug) as $debug # informative
| .;
 
# input: output of pp
# check that length is sufficient, and that $maxm and $maxn are large enough
def adequate($maxm; $maxn):
# ( ".[$N-1].m is \(.[$N-1].m) vs $maxm=\($maxm)" | debug) as $debug |
# ( ".[$N-1].n is \(.[$N-1].n) vs $maxn=\($maxn)" | debug) as $debug |
length > $N
and .[$N-1].m < $maxm - 3 # -2 is not sufficient
and .[$N-1].n < $maxn - 3 # -2 is not sufficient
;
 
# If our search has not been `adequate` then increase $maxm and $maxn
# input: [maxm, maxn, ([maxm,maxn]|pp)]
# output: pp
def adapt:
. as [$maxm, $maxn, $one]
| if ($one|adequate($maxm; $maxn)) then $one
else [$maxm + 2, $maxn + 1.6] as $maxplus
# | ("retrying with \($maxplus)" | debug) as $debug
| ($maxplus|pp) as $two
| $maxplus + [$two] | adapt
end;
 
# We start by selecting m and n so that
# m*n >> $N, i.e., 0.63 * m^2 >> $N , so m >> sqrt(1.585 * $N)
# Using 7 as the constant to start with is sufficient to avoid too much rework.
((9 * $N) | sqrt) as $maxm
| (0.63 * $maxm + 1) as $maxn
| [$maxm, $maxn] as $max
| ($max | pp) as $pp
| ($max + [$pp]) | [adapt[:$N][].p] ;
 
# The stretch goal:
def stretch:
250
| "\nThe \(.)th Pierpoint prime of the first kind is \(pierponts(true)[-1])",
"\nThe \(.)th Pierpoint prime of the second kind is \(pierponts(false)[-1])" ;
 
# The primary task:
50
| "\nThe first \(.) Pierpoint primes of the first kind:", (pierponts(true) | table(10;8)),
"\nThe first \(.) Pierpoint primes of the second kind:", (pierponts(false) | table(10;8))</syntaxhighlight>
{{out}}
<pre>
The first 50 Pierpoint primes of the first kind:
2 3 5 7 13 17 19 37 73 97
109 163 193 257 433 487 577 769 1153 1297
1459 2593 2917 3457 3889 10369 12289 17497 18433 39367
52489 65537 139969 147457 209953 331777 472393 629857 746497 786433
839809 995329 1179649 1492993 1769473 1990657 2654209 5038849 5308417 8503057
 
The first 50 Pierpoint primes of the second kind:
2 3 5 7 11 17 23 31 47 53
71 107 127 191 383 431 647 863 971 1151
2591 4373 6143 6911 8191 8747 13121 15551 23327 27647
62207 73727 131071 139967 165887 294911 314927 442367 472391 497663
524287 786431 995327 1062881 2519423 10616831 17915903 18874367 25509167 30233087
</pre>
 
=={{header|Julia}}==
The generator method is very fast but does not guarantee the primes are generated in order. Therefore we generate two times the primes needed and then sort and return the lower half.
<langsyntaxhighlight lang="julia">using Primes
 
function pierponts(N, firstkind = true)
Line 623 ⟶ 1,572:
 
println("\nThe 2000th Pierpont prime (second kind) is: ", pierponts(2000, false)[2000])
</langsyntaxhighlight>{{out}}
<pre>
The first 50 Pierpont primes (first kind) are: BigInt[2, 3, 5, 7, 13, 17, 19, 37, 73, 97, 109, 163, 193, 257, 433, 487, 577, 769, 1153, 1297, 1459, 2593, 2917, 3457, 3889, 10369, 12289, 17497, 18433, 39367, 52489, 65537, 139969, 147457, 209953, 331777, 472393, 629857, 746497, 786433, 839809, 995329, 1179649, 1492993, 1769473, 1990657, 2654209, 5038849, 5308417, 8503057]
Line 644 ⟶ 1,593:
=={{header|Kotlin}}==
{{trans|Go}}
<langsyntaxhighlight lang="scala">import java.math.BigInteger
import kotlin.math.min
 
Line 721 ⟶ 1,670:
println("\n2000th Pierpont prime of the first kind: ${p[0][1999]}")
println("\n2000th Pierpont prime of the first kind: ${p[1][1999]}")
}</langsyntaxhighlight>
{{out}}
<pre>First 50 Pierpont primes of the first kind:
Line 748 ⟶ 1,697:
 
2000th Pierpont prime of the first kind: 1702224134662426018061116932011222570937093650174807121918750428723338890211147039320296240754205680537318845776107057915956535566573559841027244444877454493022783449689509569107393738917120492483994302725479684822283929715327187974256253064796234576415398735760543848603844607</pre>
 
=={{header|Lua}}==
<syntaxhighlight lang="lua">local function isprime(n)
if n < 2 then return false end
if n % 2 == 0 then return n==2 end
if n % 3 == 0 then return n==3 end
local f, limit = 5, math.sqrt(n)
for f = 5, limit, 6 do
if n % f == 0 then return false end
if n % (f+2) == 0 then return false end
end
return true
end
 
local function s3iter()
local s, i2, i3 = {1}, 1, 1
return function()
local n2, n3, val = 2*s[i2], 3*s[i3], s[#s]
s[#s+1] = math.min(n2, n3)
i2, i3 = i2 + (n2<=n3 and 1 or 0), i3 + (n2>=n3 and 1 or 0)
return val
end
end
 
local function pierpont(n)
local list1, list2, s3next = {}, {}, s3iter()
while #list1 < n or #list2 < n do
local s3 = s3next()
if #list1 < n then
if isprime(s3+1) then list1[#list1+1] = s3+1 end
end
if #list2 < n then
if isprime(s3-1) then list2[#list2+1] = s3-1 end
end
end
return list1, list2
end
 
local N = 50
local p1, p2 = pierpont(N)
 
print("First 50 Pierpont primes of the first kind:")
for i, p in ipairs(p1) do
io.write(string.format("%8d%s", p, (i%10==0 and "\n" or " ")))
end
 
print()
 
print("First 50 Pierpont primes of the second kind:")
for i, p in ipairs(p2) do
io.write(string.format("%8d%s", p, (i%10==0 and "\n" or " ")))
end</syntaxhighlight>
{{out}}
<pre>First 50 Pierpont primes of the first kind:
2 3 5 7 13 17 19 37 73 97
109 163 193 257 433 487 577 769 1153 1297
1459 2593 2917 3457 3889 10369 12289 17497 18433 39367
52489 65537 139969 147457 209953 331777 472393 629857 746497 786433
839809 995329 1179649 1492993 1769473 1990657 2654209 5038849 5308417 8503057
 
First 50 Pierpont primes of the second kind:
2 3 5 7 11 17 23 31 47 53
71 107 127 191 383 431 647 863 971 1151
2591 4373 6143 6911 8191 8747 13121 15551 23327 27647
62207 73727 131071 139967 165887 294911 314927 442367 472391 497663
524287 786431 995327 1062881 2519423 10616831 17915903 18874367 25509167 30233087</pre>
 
=={{header|M2000 Interpreter}}==
{{trans|freebasic}}
<syntaxhighlight lang="m2000 interpreter">Module Pierpoint_Primes {
Form 80
Set Fast !
const NPP=50
dim pier(1 to 2, 1 to NPP), np(1 to 2) = 0
def long x = 1, j
while np(1)<=NPP or np(2)<=NPP
x++
j = @is_pierpont(x)
if j>0 Else Continue
if j mod 2 = 1 then np(1)++ :if np(1) <= NPP then pier(1, np(1)) = x
if j > 1 then np(2)++ : if np(2) <= NPP then pier(2, np(2)) = x
end while
print "First ";NPP;" Pierpont primes of the first kind:"
for j = 1 to NPP
print pier(2, j),
next j
if pos>0 then print
print "First ";NPP;" Pierpont primes of the second kind:"
for j = 1 to NPP
print pier(1, j),
next j
if pos>0 then print
Set Fast
function is_prime(n as decimal)
if n < 2 then = false : exit function
if n <4 then = true : exit function
if n mod 2 = 0 then = false : exit function
local i as long
for i = 3 to int(sqrt(n))+1 step 2
if n mod i = 0 then = false : exit function
next i
= true
end function
function is_23(n as long)
while n mod 2 = 0
n = n div 2
end while
while n mod 3 = 0
n = n div 3
end while
if n = 1 then = true else = false
end function
function is_pierpont(n as long)
if not @is_prime(n) then = 0& : exit function 'not prime
Local p1 = @is_23(n+1), p2 = @is_23(n-1)
if p1 and p2 then = 3 : exit function 'pierpont prime of both kinds
if p1 then = 1 : exit function 'pierpont prime of the 1st kind
if p2 then = 2 : exit function 'pierpont prime of the 2nd kind
= 0 'prime, but not pierpont
end function
}
Pierpoint_Primes</syntaxhighlight>
 
{{out}}
<pre>
First 50 Pierpont primes of the first kind:
2 3 5 7 13 17 19 37 73 97
109 163 193 257 433 487 577 769 1153 1297
1459 2593 2917 3457 3889 10369 12289 17497 18433 39367
52489 65537 139969 147457 209953 331777 472393 629857 746497 786433
839809 995329 1179649 1492993 1769473 1990657 2654209 5038849 5308417 8503057
 
First 50 Pierpont primes of the second kind:
2 3 5 7 11 17 23 31 47 53
71 107 127 191 383 431 647 863 971 1151
2591 4373 6143 6911 8191 8747 13121 15551 23327 27647
62207 73727 131071 139967 165887 294911 314927 442367 472391 497663
524287 786431 995327 1062881 2519423 10616831 17915903 18874367 25509167 30233087
</pre>
 
=={{header|Mathematica}}/{{header|Wolfram Language}}==
<syntaxhighlight lang="mathematica">ClearAll[FindUpToMax]
FindUpToMax[max_Integer, b_Integer] := Module[{res, num},
res = {};
Do[
num = 2^u 3^v + b;
If[PrimeQ[num], AppendTo[res, num]]
,
{u, 0, Ceiling@Log[2, max]}
,
{v, 0, Ceiling@Log[3, max]}
];
res //= Select[LessEqualThan[max]];
res //= Sort;
res
]
Print["Piermont primes of the first kind:"]
Take[FindUpToMax[10^10, +1], UpTo[50]]
Print["Piermont primes of the second kind:"]
Take[FindUpToMax[10^10, -1], UpTo[50]]
Print["250th Piermont prime of the first kind:"]
Part[FindUpToMax[10^34, +1], 250]
Print["250th Piermont prime of the second kind:"]
Part[FindUpToMax[10^34, -1], 250]
Print["1000th Piermont prime of the first kind:"]
Part[FindUpToMax[10^130, +1], 1000]
Print["1000th Piermont prime of the second kind:"]
Part[FindUpToMax[10^150, -1], 1000]</syntaxhighlight>
{{out}}
<pre>Piermont primes of the first kind:
{2,3,5,7,13,17,19,37,73,97,109,163,193,257,433,487,577,769,1153,1297,1459,2593,2917,3457,3889,10369,12289,17497,18433,39367,52489,65537,139969,147457,209953,331777,472393,629857,746497,786433,839809,995329,1179649,1492993,1769473,1990657,2654209,5038849,5308417,8503057}
Piermont primes of the second kind:
{2,3,5,7,11,17,23,31,47,53,71,107,127,191,383,431,647,863,971,1151,2591,4373,6143,6911,8191,8747,13121,15551,23327,27647,62207,73727,131071,139967,165887,294911,314927,442367,472391,497663,524287,786431,995327,1062881,2519423,10616831,17915903,18874367,25509167,30233087}
250th Piermont primes of the first kind:
62518864539857068333550694039553
250th Piermont primes of the second kind:
4111131172000956525894875083702271
1000th Piermont prime of the first kind:
69269314716439690250482558089997110961545818230232043107188537422260188701607997086273960899938499201024414931399264696270849
1000th Piermont prime of the second kind:
1308088756227965581249669045506775407896673213729433892383353027814827286537163695213418982500477392209371001259166465228280492460735463423</pre>
 
 
=={{header|Nim}}==
<syntaxhighlight lang="nim">import math, strutils
 
func isPrime(n: int): bool =
## Check if "n" is prime by trying successive divisions.
## "n" is supposed not to be a multiple of 2 or 3.
var d = 5
var delta = 2
while d <= int(sqrt(n.toFloat)):
if n mod d == 0: return false
inc d, delta
delta = 6 - delta
result = true
 
func isProduct23(n: int): bool =
## Check if "n" has only 2 and 3 for prime divisors
## (i.e. that "n = 2^u * 3^v").
var n = n
while (n and 1) == 0: n = n shr 1
while n mod 3 == 0: n = n div 3
result = n == 1
 
iterator pierpont(maxCount: Positive; k: int): int =
## Yield the Pierpoint primes of first or second kind according
## to the value of "k" (+1 for first kind, -1 for second kind).
yield 2
yield 3
var n = 5
var delta = 2 # 2 and 4 alternatively to skip the multiples of 2 and 3.
yield n
var count = 3
while count < maxCount:
inc n, delta
delta = 6 - delta
if isProduct23(n - k) and isPrime(n):
inc count
yield n
 
template pierpont1*(maxCount: Positive): int = pierpont(maxCount, +1)
template pierpont2*(maxCount: Positive): int = pierpont(maxCount, -1)
 
 
when isMainModule:
 
echo "First 50 Pierpont primes of the first kind:"
var count = 0
for n in pierpont1(50):
stdout.write ($n).align(9)
inc count
if count mod 10 == 0: stdout.write '\n'
 
echo ""
echo "First 50 Pierpont primes of the second kind:"
count = 0
for n in pierpont2(50):
stdout.write ($n).align(9)
inc count
if count mod 10 == 0: stdout.write '\n'</syntaxhighlight>
 
{{out}}
<pre>First 50 Pierpont primes of the first kind:
2 3 5 7 13 17 19 37 73 97
109 163 193 257 433 487 577 769 1153 1297
1459 2593 2917 3457 3889 10369 12289 17497 18433 39367
52489 65537 139969 147457 209953 331777 472393 629857 746497 786433
839809 995329 1179649 1492993 1769473 1990657 2654209 5038849 5308417 8503057
 
First 50 Pierpont primes of the second kind:
2 3 5 7 11 17 23 31 47 53
71 107 127 191 383 431 647 863 971 1151
2591 4373 6143 6911 8191 8747 13121 15551 23327 27647
62207 73727 131071 139967 165887 294911 314927 442367 472391 497663
524287 786431 995327 1062881 2519423 10616831 17915903 18874367 25509167 30233087</pre>
 
=={{header|Perl}}==
{{trans|Raku}}
{{libheader|ntheory}}
<langsyntaxhighlight lang="perl">use strict;
use warnings;
use feature 'say';
Line 804 ⟶ 2,012:
 
say "\n250th Pierpont prime of the first kind: " . $pierpont_1st[249];
say "\n250th Pierpont prime of the second kind: " . $pierpont_2nd[249];</langsyntaxhighlight>
{{out}}
<pre>First 50 Pierpont primes of the first kind:
Line 825 ⟶ 2,033:
 
=={{header|Phix}}==
{{libheader|Phix/mpfr}}
{{trans|Go}}
I also tried using a priority queue, popping the smallest (and any duplicates of it) and pushing *2 and *3 on every iteration, which worked pretty well but ended up about 20% slower, with about 1500 untested candidates left in the queue by the time it found the 2000th second kind.
<!--<syntaxhighlight lang="phix">(phixonline)-->
<lang Phix>-- demo/rosetta/Pierpont_primes.exw
<span style="color: #000080;font-style:italic;">-- demo/rosetta/Pierpont_primes.exw</span>
include mpfr.e
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span>
<span style="color: #008080;">include</span> <span style="color: #004080;">mpfr</span><span style="color: #0000FF;">.</span><span style="color: #000000;">e</span>
function pierpont(integer n)
sequence p = {{mpz_init(2)},{}},
<span style="color: #008080;">function</span> <span style="color: #000000;">pierpont</span><span style="color: #0000FF;">(</span><span style="color: #004080;">integer</span> <span style="color: #000000;">n</span><span style="color: #0000FF;">)</span>
s = {mpz_init(1)}
<span style="color: #004080;">sequence</span> <span style="color: #000000;">p</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{{</span><span style="color: #7060A8;">mpz_init</span><span style="color: #0000FF;">(</span><span style="color: #000000;">2</span><span style="color: #0000FF;">)},{}},</span>
integer i2 = 1, i3 = 1
<span style="color: #000000;">s</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{</span><span style="color: #7060A8;">mpz_init</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">)}</span>
mpz {n2, n3, t} = mpz_inits(3)
<span style="color: #004080;">integer</span> <span style="color: #000000;">i2</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">i3</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">1</span>
randstate state = gmp_randinit_mt()
<span style="color: #004080;">mpz</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">n2</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">n3</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">t</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">mpz_inits</span><span style="color: #0000FF;">(</span><span style="color: #000000;">3</span><span style="color: #0000FF;">)</span>
atom t1 = time()+1
<span style="color: #004080;">atom</span> <span style="color: #000000;">t1</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">time</span><span style="color: #0000FF;">()+</span><span style="color: #000000;">1</span>
while min(length(p[1]),length(p[2])) < n do
<span style="color: #008080;">while</span> <span style="color: #7060A8;">min</span><span style="color: #0000FF;">(</span><span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">p</span><span style="color: #0000FF;">[</span><span style="color: #000000;">1</span><span style="color: #0000FF;">]),</span><span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">p</span><span style="color: #0000FF;">[</span><span style="color: #000000;">2</span><span style="color: #0000FF;">]))</span> <span style="color: #0000FF;"><</span> <span style="color: #000000;">n</span> <span style="color: #008080;">do</span>
mpz_mul_si(n2,s[i2],2)
<span style="color: #7060A8;">mpz_mul_si</span><span style="color: #0000FF;">(</span><span style="color: #000000;">n2</span><span style="color: #0000FF;">,</span><span style="color: #000000;">s</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i2</span><span style="color: #0000FF;">],</span><span style="color: #000000;">2</span><span style="color: #0000FF;">)</span>
mpz_mul_si(n3,s[i3],3)
<span style="color: #7060A8;">mpz_mul_si</span><span style="color: #0000FF;">(</span><span style="color: #000000;">n3</span><span style="color: #0000FF;">,</span><span style="color: #000000;">s</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i3</span><span style="color: #0000FF;">],</span><span style="color: #000000;">3</span><span style="color: #0000FF;">)</span>
if mpz_cmp(n2,n3)<0 then
<span style="color: #008080;">if</span> <span style="color: #7060A8;">mpz_cmp</span><span style="color: #0000FF;">(</span><span style="color: #000000;">n2</span><span style="color: #0000FF;">,</span><span style="color: #000000;">n3</span><span style="color: #0000FF;">)<</span><span style="color: #000000;">0</span> <span style="color: #008080;">then</span>
mpz_set(t,n2)
<span style="color: #7060A8;">mpz_set</span><span style="color: #0000FF;">(</span><span style="color: #000000;">t</span><span style="color: #0000FF;">,</span><span style="color: #000000;">n2</span><span style="color: #0000FF;">)</span>
i2 += 1
<span style="color: #000000;">i2</span> <span style="color: #0000FF;">+=</span> <span style="color: #000000;">1</span>
else
<span style="color: mpz_set(t,n3)#008080;">else</span>
<span style="color: #7060A8;">mpz_set</span><span style="color: #0000FF;">(</span><span style="color: #000000;">t</span><span style="color: #0000FF;">,</span><span style="color: #000000;">n3</span><span style="color: #0000FF;">)</span>
i3 += 1
<span style="color: #000000;">i3</span> <span style="color: #0000FF;">+=</span> <span style="color: #000000;">1</span>
end if
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
if mpz_cmp(t,s[$]) > 0 then
<span style="color: #008080;">if</span> <span style="color: #7060A8;">mpz_cmp</span><span style="color: #0000FF;">(</span><span style="color: #000000;">t</span><span style="color: #0000FF;">,</span><span style="color: #000000;">s</span><span style="color: #0000FF;">[$])</span> <span style="color: #0000FF;">></span> <span style="color: #000000;">0</span> <span style="color: #008080;">then</span>
s = append(s, mpz_init_set(t))
<span style="color: #000000;">s</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">append</span><span style="color: #0000FF;">(</span><span style="color: #000000;">s</span><span style="color: #0000FF;">,</span> <span style="color: #7060A8;">mpz_init_set</span><span style="color: #0000FF;">(</span><span style="color: #000000;">t</span><span style="color: #0000FF;">))</span>
mpz_add_ui(t,t,1)
<span style="color: #7060A8;">mpz_add_ui</span><span style="color: #0000FF;">(</span><span style="color: #000000;">t</span><span style="color: #0000FF;">,</span><span style="color: #000000;">t</span><span style="color: #0000FF;">,</span><span style="color: #000000;">1</span><span style="color: #0000FF;">)</span>
if length(p[1]) < n and mpz_probable_prime_p(t,state) then
<span style="color: #008080;">if</span> <span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">p</span><span style="color: #0000FF;">[</span><span style="color: #000000;">1</span><span style="color: #0000FF;">])</span> <span style="color: #0000FF;"><</span> <span style="color: #000000;">n</span> <span style="color: #008080;">and</span> <span style="color: #7060A8;">mpz_prime</span><span style="color: #0000FF;">(</span><span style="color: #000000;">t</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">then</span>
p[1] = append(p[1],mpz_init_set(t))
<span style="color: #000000;">p</span><span style="color: #0000FF;">[</span><span style="color: #000000;">1</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">append</span><span style="color: #0000FF;">(</span><span style="color: #000000;">p</span><span style="color: #0000FF;">[</span><span style="color: #000000;">1</span><span style="color: #0000FF;">],</span><span style="color: #7060A8;">mpz_init_set</span><span style="color: #0000FF;">(</span><span style="color: #000000;">t</span><span style="color: #0000FF;">))</span>
end if
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
mpz_sub_ui(t,t,2)
<span style="color: #7060A8;">mpz_sub_ui</span><span style="color: #0000FF;">(</span><span style="color: #000000;">t</span><span style="color: #0000FF;">,</span><span style="color: #000000;">t</span><span style="color: #0000FF;">,</span><span style="color: #000000;">2</span><span style="color: #0000FF;">)</span>
if length(p[2]) < n and mpz_probable_prime_p(t,state) then
<span style="color: #008080;">if</span> <span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">p</span><span style="color: #0000FF;">[</span><span style="color: #000000;">2</span><span style="color: #0000FF;">])</span> <span style="color: #0000FF;"><</span> <span style="color: #000000;">n</span> <span style="color: #008080;">and</span> <span style="color: #7060A8;">mpz_prime</span><span style="color: #0000FF;">(</span><span style="color: #000000;">t</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">then</span>
p[2] = append(p[2],mpz_init_set(t))
<span style="color: #000000;">p</span><span style="color: #0000FF;">[</span><span style="color: #000000;">2</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">append</span><span style="color: #0000FF;">(</span><span style="color: #000000;">p</span><span style="color: #0000FF;">[</span><span style="color: #000000;">2</span><span style="color: #0000FF;">],</span><span style="color: #7060A8;">mpz_init_set</span><span style="color: #0000FF;">(</span><span style="color: #000000;">t</span><span style="color: #0000FF;">))</span>
end if
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
end if
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
if time()>t1 then
<span style="color: #008080;">if</span> <span style="color: #7060A8;">time</span><span style="color: #0000FF;">()></span><span style="color: #000000;">t1</span> <span style="color: #008080;">then</span>
printf(1,"Found: 1st:%d/%d, 2nd:%d/%d\r",
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"Found: 1st:%d/%d, 2nd:%d/%d\r"</span><span style="color: #0000FF;">,</span>
{length(p[1]),n,length(p[2]),n})
<span style="color: #0000FF;">{</span><span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">p</span><span style="color: #0000FF;">[</span><span style="color: #000000;">1</span><span style="color: #0000FF;">]),</span><span style="color: #000000;">n</span><span style="color: #0000FF;">,</span><span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">p</span><span style="color: #0000FF;">[</span><span style="color: #000000;">2</span><span style="color: #0000FF;">]),</span><span style="color: #000000;">n</span><span style="color: #0000FF;">})</span>
t1 = time()+1
<span style="color: #000000;">t1</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">time</span><span style="color: #0000FF;">()+</span><span style="color: #000000;">1</span>
end if
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
end while
<span style="color: #008080;">end</span> <span style="color: #008080;">while</span>
return p
<span style="color: #008080;">return</span> <span style="color: #000000;">p</span>
end function
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
constant limit = 2000 -- 2 mins
<span style="color: #000080;font-style:italic;">--constant limit = 10002000 -- 8.1s2 mins
--constant limit = 250 1000 -- 08.1s</span>
<span style="color: #008080;">constant</span> <span style="color: #000000;">limit</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">250</span> <span style="color: #000080;font-style:italic;">-- 0.1s</span>
atom t0 = time()
<span style="color: #004080;">atom</span> <span style="color: #000000;">t0</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">time</span><span style="color: #0000FF;">()</span>
sequence p = pierpont(limit)
<span style="color: #004080;">sequence</span> <span style="color: #000000;">p</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">pierpont</span><span style="color: #0000FF;">(</span><span style="color: #000000;">limit</span><span style="color: #0000FF;">)</span>
constant fs = {"first","second"}
<span style="color: #008080;">constant</span> <span style="color: #000000;">fs</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{</span><span style="color: #008000;">"first"</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"second"</span><span style="color: #0000FF;">}</span>
 
for i=1 to length(fs) do
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">fs</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">do</span>
printf(1,"First 50 Pierpont primes of the %s kind:\n",{fs[i]})
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"First 50 Pierpont primes of the %s kind:\n"</span><span style="color: #0000FF;">,{</span><span style="color: #000000;">fs</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">]})</span>
for j=1 to 50 do
<span style="color: #008080;">for</span> <span style="color: #000000;">j</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #000000;">50</span> <span style="color: #008080;">do</span>
mpfr_printf(1,"%8Zd ", p[i][j])
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"%8s "</span><span style="color: #0000FF;">,</span> <span style="color: #0000FF;">{</span><span style="color: #7060A8;">mpz_get_str</span><span style="color: #0000FF;">(</span><span style="color: #000000;">p</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">][</span><span style="color: #000000;">j</span><span style="color: #0000FF;">])})</span>
if mod(j,10)=0 then printf(1,"\n") end if
<span style="color: #008080;">if</span> <span style="color: #7060A8;">mod</span><span style="color: #0000FF;">(</span><span style="color: #000000;">j</span><span style="color: #0000FF;">,</span><span style="color: #000000;">10</span><span style="color: #0000FF;">)=</span><span style="color: #000000;">0</span> <span style="color: #008080;">then</span> <span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"\n"</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
end for
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
printf(1,"\n")
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"\n"</span><span style="color: #0000FF;">)</span>
end for
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
constant t = {250,1000,2000}
<span style="color: #008080;">constant</span> <span style="color: #000000;">t</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">250</span><span style="color: #0000FF;">,</span><span style="color: #000000;">1000</span><span style="color: #0000FF;">,</span><span style="color: #000000;">2000</span><span style="color: #0000FF;">}</span>
for i=1 to length(t) do
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">t</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">do</span>
integer ti = t[i]
<span style="color: #004080;">integer</span> <span style="color: #000000;">ti</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">t</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">]</span>
if ti>limit then exit end if
<span style="color: #008080;">if</span> <span style="color: #000000;">ti</span><span style="color: #0000FF;">></span><span style="color: #000000;">limit</span> <span style="color: #008080;">then</span> <span style="color: #008080;">exit</span> <span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
for j=1 to length(fs) do
<span style="color: #008080;">for</span> <span style="color: #000000;">j</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">fs</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">do</span>
string zs = shorten(mpz_get_str(p[j][ti]))
<span style="color: #004080;">string</span> <span style="color: #000000;">zs</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">shorten</span><span style="color: #0000FF;">(</span><span style="color: #7060A8;">mpz_get_str</span><span style="color: #0000FF;">(</span><span style="color: #000000;">p</span><span style="color: #0000FF;">[</span><span style="color: #000000;">j</span><span style="color: #0000FF;">][</span><span style="color: #000000;">ti</span><span style="color: #0000FF;">]))</span>
printf(1,"%dth Pierpont prime of the %s kind: %s\n",{ti,fs[j],zs})
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"%dth Pierpont prime of the %s kind: %s\n"</span><span style="color: #0000FF;">,{</span><span style="color: #000000;">ti</span><span style="color: #0000FF;">,</span><span style="color: #000000;">fs</span><span style="color: #0000FF;">[</span><span style="color: #000000;">j</span><span style="color: #0000FF;">],</span><span style="color: #000000;">zs</span><span style="color: #0000FF;">})</span>
end for
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
printf(1,"\n")
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"\n"</span><span style="color: #0000FF;">)</span>
end for
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
printf(1,"Took %s\n", elapsed(time()-t0))</lang>
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"Took %s\n"</span><span style="color: #0000FF;">,</span> <span style="color: #7060A8;">elapsed</span><span style="color: #0000FF;">(</span><span style="color: #7060A8;">time</span><span style="color: #0000FF;">()-</span><span style="color: #000000;">t0</span><span style="color: #0000FF;">))</span>
<!--</syntaxhighlight>-->
{{out}}
<pre>
Line 922 ⟶ 2,132:
Took 2 minutes and 01s
</pre>
Note that shorten() has recently been added as a builtin, in honour of this and the several other dozen rc
tasks that previously all re-implemented some variation of it. Should you not yet have it, just use this:
<lang Phix>function shorten(string s, what="digits", integer ml=20)
integer l = length(s)
string ls = sprintf(" (%d %s)",{l,what})
if l>ml*2+3+length(ls) then
s[ml..-ml] = "..."
s &= ls
end if
return s
end function</lang>
 
=={{header|Prolog}}==
<syntaxhighlight lang="prolog">
?- use_module(library(heaps)).
 
three_smooth(Lz) :-
singleton_heap(H, 1, nothing),
lazy_list(next_3smooth, 0-H, Lz).
 
next_3smooth(Top-H0, N-H2, N) :-
min_of_heap(H0, Top, _), !,
get_from_heap(H0, Top, _, H1),
next_3smooth(Top-H1, N-H2, N).
next_3smooth(_-H0, N-H3, N) :-
get_from_heap(H0, N, _, H1),
N2 is N * 2,
N3 is N * 3,
add_to_heap(H1, N2, nothing, H2),
add_to_heap(H2, N3, nothing, H3).
 
first_kind(K) :-
three_smooth(Ns), member(N, Ns),
K is N + 1,
prime(K).
 
second_kind(K) :-
three_smooth(Ns), member(N, Ns),
K is N - 1,
prime(K).
 
show(Seq, N) :-
format("The first ~w values of ~s are: ", [N, Seq]),
once(findnsols(N, X, call(Seq, X), L)),
write(L), nl,
once(offset(249, call(Seq, TwoFifty))),
format("The 250th value of ~w is ~w~n", [Seq, TwoFifty]).
 
main :-
show(first_kind, 50), nl,
show(second_kind, 50), nl,
halt.
 
% primality checker -- Miller Rabin preceded with a round of trial divisions.
 
prime(N) :-
integer(N),
N > 1,
divcheck(
N,
[ 2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31,
37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79,
83, 89, 97, 101, 103, 107, 109, 113, 127, 131, 137,
139, 149],
Result),
((Result = prime, !); miller_rabin_primality_test(N)).
 
divcheck(_, [], unknown) :- !.
divcheck(N, [P|_], prime) :- P*P > N, !.
divcheck(N, [P|Ps], State) :- N mod P =\= 0, divcheck(N, Ps, State).
 
miller_rabin_primality_test(N) :-
bases(Bases, N),
forall(member(A, Bases), strong_fermat_pseudoprime(N, A)).
 
miller_rabin_precision(16).
 
bases([31, 73], N) :- N < 9_080_191, !.
bases([2, 7, 61], N) :- N < 4_759_123_141, !.
bases([2, 325, 9_375, 28_178, 450_775, 9_780_504, 1_795_265_022], N) :-
N < 18_446_744_073_709_551_616, !. % 2^64
bases(Bases, N) :-
miller_rabin_precision(T), RndLimit is N - 2,
length(Bases, T), maplist(random_between(2, RndLimit), Bases).
 
strong_fermat_pseudoprime(N, A) :- % miller-rabin strong pseudoprime test with base A.
succ(Pn, N), factor_2s(Pn, S, D),
X is powm(A, D, N),
((X =:= 1, !); \+ composite_witness(N, S, X)).
 
composite_witness(_, 0, _) :- !.
composite_witness(N, K, X) :-
X =\= N-1,
succ(Pk, K), X2 is (X*X) mod N, composite_witness(N, Pk, X2).
 
factor_2s(N, S, D) :- factor_2s(0, N, S, D).
factor_2s(S, D, S, D) :- D /\ 1 =\= 0, !.
factor_2s(S0, D0, S, D) :-
succ(S0, S1), D1 is D0 >> 1,
factor_2s(S1, D1, S, D).
 
?- main.
</syntaxhighlight>
{{Out}}
<pre>
The first 50 values of first_kind are: [2,3,5,7,13,17,19,37,73,97,109,163,193,257,433,487,577,769,1153,1297,1459,2593,2917,3457,3889,10369,12289,17497,18433,39367,52489,65537,139969,147457,209953,331777,472393,629857,746497,786433,839809,995329,1179649,1492993,1769473,1990657,2654209,5038849,5308417,8503057]
The 250th value of first_kind is 62518864539857068333550694039553
 
The first 50 values of second_kind are: [2,3,5,7,11,17,23,31,47,53,71,107,127,191,383,431,647,863,971,1151,2591,4373,6143,6911,8191,8747,13121,15551,23327,27647,62207,73727,131071,139967,165887,294911,314927,442367,472391,497663,524287,786431,995327,1062881,2519423,10616831,17915903,18874367,25509167,30233087]
The 250th value of second_kind is 4111131172000956525894875083702271
</pre>
=={{header|Python}}==
{{trans|Go}}
<langsyntaxhighlight lang="python">import random
 
# Copied from https://rosettacode.org/wiki/Miller-Rabin_primality_test#Python
Line 1,013 ⟶ 2,311:
print "250th Pierpont prime of the second kind:", pp2[249]
 
main()</langsyntaxhighlight>
{{out}}
<pre>First 50 Pierpont primes of the first kind:
Line 1,029 ⟶ 2,327:
250th Pierpont prime of the first kind: 62518864539857068333550694039553
250th Pierpont prime of the second kind: 4111131172000956525894875083702271</pre>
 
=={{header|Quackery}}==
 
Uses <code>smoothwith</code> from [[N-smooth numbers#Quackery]] and <code>isprime</code> from [[Primality by trial division#Quackery]].
 
Using trial division to determine primality makes finding the 250th Pierpont primes impractical. This should be updated when a coding of a more suitable test becomes available.
 
<syntaxhighlight lang="quackery"> [ stack ] is pierponts
[ stack ] is kind
[ stack ] is quantity
 
[ 1 - -2 * 1+ kind put
1+ quantity put
' [ -1 ] pierponts put
' [ 2 3 ] smoothwith
[ -1 peek
kind share +
dup isprime not iff
[ drop false ] done
pierponts share -1 peek
over = iff
[ drop false ] done
pierponts take
swap join
dup size swap
pierponts put
quantity share = ]
drop
quantity release
kind release
pierponts take
behead drop ] is pierpontprimes
 
say "Pierpont primes of the first kind." cr
50 1 pierpontprimes echo
cr cr
say "Pierpont primes of the second kind." cr
50 2 pierpontprimes echo</syntaxhighlight>
 
{{out}}
 
<pre>Pierpont primes of the first kind.
[ 2 3 5 7 13 17 19 37 73 97 109 163 193 257 433 487 577 769 1153 1297 1459 2593 2917 3457 3889 10369 12289 17497 18433 39367 52489 65537 139969 147457 209953 331777 472393 629857 746497 786433 839809 995329 1179649 1492993 1769473 1990657 2654209 5038849 5308417 8503057 ]
 
Pierpont primes of the second kind.
[ 2 3 5 7 11 17 23 31 47 53 71 107 127 191 383 431 647 863 971 1151 2591 4373 6143 6911 8191 8747 13121 15551 23327 27647 62207 73727 131071 139967 165887 294911 314927 442367 472391 497663 524287 786431 995327 1062881 2519423 10616831 17915903 18874367 25509167 30233087 ]</pre>
 
=={{header|Raku}}==
Line 1,040 ⟶ 2,384:
an overabundance. No need to rely on magic numbers. No need to sort them. It
produces exactly what is needed, in order, on demand.
{{libheader|ntheory}}
 
<syntaxhighlight lang="raku" perl6line>use ntheory:from<Perl5> <is_prime>;
 
sub pierpont ($kind is copy = 1) {
Line 1,074 ⟶ 2,418:
say "\n250th Pierpont prime of the first kind: " ~ pierpont[249];
 
say "\n250th Pierpont prime of the second kind: " ~ pierpont(2)[249];</langsyntaxhighlight>
 
{{out}}
Line 1,100 ⟶ 2,444:
(Cut down output as it is exactly the same as the first version for {2,3} +1 and {2,3} -1; leaves room to demo some other options.)
 
<syntaxhighlight lang="raku" perl6line>sub smooth-numbers (*@list) {
cache my \Smooth := gather {
my %i = (flat @list) Z=> (Smooth.iterator for ^@list);
Line 1,128 ⟶ 2,472:
"\nOEIS: A077314:", (2,7), -1, 20,
"\nOEIS: A174144 - (\"Humble\" primes 1st):", (2,3,5,7), 1, 20,
"\nOEIS: A299171A347977 - (\"Humble\" primes 2nd):", (2,3,5,7), -1, 20,
"\nOEIS: A077499:", (2,11), 1, 20,
"\nOEIS: A077315:", (2,11), -1, 20,
Line 1,138 ⟶ 2,482:
say "$title smooth \{$primes\} {$add > 0 ?? '+' !! '-'} 1 ";
put smooth-numbers(|$primes).map( * + $add ).grep( *.is-prime )[^$count]
}</langsyntaxhighlight>
 
{{Out}}
Line 1,192 ⟶ 2,536:
The REXX language has a "big num" capability to handle almost any amount of decimal digits, &nbsp; but
<br>it lacks a robust &nbsp; '''isPrime''' &nbsp; function. &nbsp; Without that, verifying very large primes is problematic.
<langsyntaxhighlight lang="rexx">/*REXX program finds and displays Pierpont primes of the first and second kinds. */
parse arg n . /*obtain optional argument from the CL.*/
if n=='' | n=="," then n= 50 /*Not specified? Then use the default.*/
numeric digits n /*ensure enough decimal digs (bit int).*/
big= copies(9, digits() ) /*BIG: used as a max number (a limit).*/
@.= '2nd'; do t=1 to -1 by -2; usum= 0; vsum= 0; s@.1= 0'1st'
#= 0do t=1 to -1 by -2; usum= 0; vsum= 0; s= 0 /*T is 1, then /*number of Pierpont primes (so far)-1. */
w #= 0 /*thenumber maxof widthPierpont ofprimes a Pierpont(so primefar). */
$=; do j=0 until #>=n /*$: the list of the" " " " Pierpont primes. */
if usum<=s then usum= get(2, 3); if vsum<=s then vsum= get(3, 2)
s= min(vsum, usum); if \isPrime(s) then iterate /*get min; isNot prime? */
#= # + 1; $= $ s /*bump counter; append.*/
w= max(w, length(s) ) end /*find max prime width.j*/
end /*j*/say
w= length(word($, #) ) /*biggest prime length.*/
say
if t==1say center(n then @= '1st'" Pierpont primes of the " @.t ' kind', max(10 /*choose(w+1), word80), for type.*/"═")
 
else @= '2nd' /* " " " " */
say center(n do "p=1 Pierpont primesby of10 the "to #; @ ' kind',_=; do k=p for max(10*; _= _ right(w+1) -1word($, 79k), "═" w)
call show $ end /*display the primes. k*/
if _\=='' then say substr( strip(_, "T"), 2)
end /*type*/
end /*p*/
exit /*stick a fork in it, we're all done. */
end /*t*/
/*──────────────────────────────────────────────────────────────────────────────────────*/
exit 0 /*stick a fork in it, we're all done. */
show: do j=1 by 10 to words($); _=
do k=j for 10; _= _ right( word($, k), w)
end /*k*/
if _\=='' then say substr( strip(_, 'T'), 2)
end /*j*/; return
/*──────────────────────────────────────────────────────────────────────────────────────*/
isPrime: procedure; parse arg x '' -1 _; if x<217 then return 0wordpos(x,"2 3 5 7 11 /*not a prime.*/13")>0
if wordpos(x,_==5 '2 3then 5return 7')\==0; then return 1 if x//2==0 then return 0 /*it's anot prime. */
if x//23==0 then return 0; if x//37==0 then return 0 /*not " " a prime.*/
do j=511 by 6 until j*j>x /*skip ÷ 3's.*/
if x//j==0 then return 0; if x//(j+2)==0 then return 0 /*not a prime. */
end /*j*/; return 1 return 1 /*it's a prime.*/
/*──────────────────────────────────────────────────────────────────────────────────────*/
get: parse arg c1,c2; m=big; do ju=0; pu= c1**ju; if pu+t>s then return min(m, pu+t)
do jv=0; pv= c2**jv; if pv >s then iterate ju
_= pu*pv + t; if _ >s then m= min(_, m)
end /*jv*/
end /*ju*/ /*see the RETURN (above). */</langsyntaxhighlight>
{{out|output|text=&nbsp; when using the default input:}}
<pre>
═════════════════════50 Pierpont primes of the 1st kind═════════════════════kind══════════════════════
2 3 5 7 13 17 19 37 73 97
109 163 193 257 433 487 577 769 1153 1297
Line 1,241 ⟶ 2,581:
839809 995329 1179649 1492993 1769473 1990657 2654209 5038849 5308417 8503057
 
══════════════════════════50 Pierpont primes of the 2nd kind══════════════════════════kind═══════════════════════════
2 3 5 7 11 17 23 31 47 53
71 107 127 191 383 431 647 863 971 1151
Line 1,249 ⟶ 2,589:
</pre>
 
=={{header|Ring}}==
<syntaxhighlight lang="ring">
load "stdlib.ring"
 
see "working..." + nl
 
pierpont = []
limit1 = 18
limit2 = 8505000
limit3 = 50
limit4 = 21
limit5 = 30500000
 
for n = 0 to limit1
for m = 0 to limit1
num = pow(2,n)*pow(3,m) + 1
if isprime(num) and num < limit2
add(pierpont,num)
ok
if num > limit2
exit
ok
next
next
 
pierpont = sort(pierpont)
 
see "First 50 Pierpont primes of the first kind:" + nl + nl
 
for n = 1 to limit3
see "" + n + ". " + pierpont[n] + nl
next
 
see "done1..." + nl
 
pierpont = []
 
for n = 0 to limit4
for m = 0 to limit4
num = pow(2,n)*pow(3,m) - 1
if isprime(num) and num < limit5
add(pierpont,num)
ok
if num > limit5
exit
ok
next
next
 
pierpont = sort(pierpont)
 
see "First 50 Pierpont primes of the second kind:" + nl + nl
 
for n = 1 to limit3
see "" + n + ". " + pierpont[n] + nl
next
 
see "done2..." + nl
</syntaxhighlight>
Output:
<pre>
working...
First 50 Pierpont primes of the first kind:
 
1. 2
2. 3
3. 5
4. 7
5. 13
6. 17
7. 19
8. 37
9. 73
10. 97
11. 109
12. 163
13. 193
14. 257
15. 433
16. 487
17. 577
18. 769
19. 1153
20. 1297
21. 1459
22. 2593
23. 2917
24. 3457
25. 3889
26. 10369
27. 12289
28. 17497
29. 18433
30. 39367
31. 52489
32. 65537
33. 139969
34. 147457
35. 209953
36. 331777
37. 472393
38. 629857
39. 746497
40. 786433
41. 839809
42. 995329
43. 1179649
44. 1492993
45. 1769473
46. 1990657
47. 2654209
48. 5038849
49. 5308417
50. 8503057
done1...
First 50 Pierpont primes of the second kind:
 
1. 2
2. 3
3. 5
4. 7
5. 11
6. 17
7. 23
8. 31
9. 47
10. 53
11. 71
12. 107
13. 127
14. 191
15. 383
16. 431
17. 647
18. 863
19. 971
20. 1151
21. 2591
22. 4373
23. 6143
24. 6911
25. 8191
26. 8747
27. 13121
28. 15551
29. 23327
30. 27647
31. 62207
32. 73727
33. 131071
34. 139967
35. 165887
36. 294911
37. 314927
38. 442367
39. 472391
40. 497663
41. 524287
42. 786431
43. 995327
44. 1062881
45. 2519423
46. 10616831
47. 17915903
48. 18874367
49. 25509167
50. 30233087
done2...
</pre>
 
=={{header|Ruby}}==
<syntaxhighlight lang="ruby">require 'gmp'
 
def smooth_generator(ar)
return to_enum(__method__, ar) unless block_given?
next_smooth = 1
queues = ar.map{|num| [num, []] }
loop do
yield next_smooth
queues.each {|m, queue| queue << next_smooth * m}
next_smooth = queues.collect{|m, queue| queue.first}.min
queues.each{|m, queue| queue.shift if queue.first == next_smooth }
end
end
def pierpont(num = 1)
return to_enum(__method__, num) unless block_given?
smooth_generator([2,3]).each{|smooth| yield smooth+num if GMP::Z(smooth + num).probab_prime? > 0}
end
 
def puts_cols(ar, n=10)
ar.each_slice(n).map{|slice|puts slice.map{|n| n.to_s.rjust(10)}.join }
end
 
n, m = 50, 250
puts "First #{n} Pierpont primes of the first kind:"
puts_cols(pierpont.take(n))
puts "#{m}th Pierpont prime of the first kind: #{pierpont.take(250).last}",""
puts "First #{n} Pierpont primes of the second kind:"
puts_cols(pierpont(-1).take(n))
puts "#{m}th Pierpont prime of the second kind: #{pierpont(-1).take(250).last}"
</syntaxhighlight>
{{out}}
<pre>First 50 Pierpont primes of the first kind:
2 3 5 7 13 17 19 37 73 97
109 163 193 257 433 487 577 769 1153 1297
1459 2593 2917 3457 3889 10369 12289 17497 18433 39367
52489 65537 139969 147457 209953 331777 472393 629857 746497 786433
839809 995329 1179649 1492993 1769473 1990657 2654209 5038849 5308417 8503057
250th Pierpont prime of the first kind: 62518864539857068333550694039553
 
First 50 Pierpont primes of the second kind:
2 3 5 7 11 17 23 31 47 53
71 107 127 191 383 431 647 863 971 1151
2591 4373 6143 6911 8191 8747 13121 15551 23327 27647
62207 73727 131071 139967 165887 294911 314927 442367 472391 497663
524287 786431 995327 1062881 2519423 10616831 17915903 18874367 25509167 30233087
250th Pierpont prime of the second kind: 4111131172000956525894875083702271
</pre>
=={{header|Sidef}}==
<langsyntaxhighlight lang="ruby">func smooth_generator(primes) {
var s = primes.len.of { [1] }
{
Line 1,278 ⟶ 2,837:
say "\n#{n}th Pierpont prime of the 1st kind: #{p}"
say "#{n}th Pierpont prime of the 2nd kind: #{q}"
}</langsyntaxhighlight>
{{out}}
<pre>
Line 1,295 ⟶ 2,854:
1000th Pierpont prime of the 1st kind: 69269314716439690250482558089997110961545818230232043107188537422260188701607997086273960899938499201024414931399264696270849
1000th Pierpont prime of the 2nd kind: 1308088756227965581249669045506775407896673213729433892383353027814827286537163695213418982500477392209371001259166465228280492460735463423
</pre>
 
=={{header|Swift}}==
 
3-smooth version.
 
Using AttaSwift's BigInt library.
 
<syntaxhighlight lang="swift">import BigInt
import Foundation
 
public func pierpoint(n: Int) -> (first: [BigInt], second: [BigInt]) {
var primes = (first: [BigInt](repeating: 0, count: n), second: [BigInt](repeating: 0, count: n))
 
primes.first[0] = 2
 
var count1 = 1, count2 = 0
var s = [BigInt(1)]
var i2 = 0, i3 = 0, k = 1
var n2 = BigInt(0), n3 = BigInt(0), t = BigInt(0)
 
while min(count1, count2) < n {
n2 = s[i2] * 2
n3 = s[i3] * 3
 
if n2 < n3 {
t = n2
i2 += 1
} else {
t = n3
i3 += 1
}
 
if t <= s[k - 1] {
continue
}
 
s.append(t)
k += 1
t += 1
 
if count1 < n && t.isPrime(rounds: 10) {
primes.first[count1] = t
count1 += 1
}
 
t -= 2
 
if count2 < n && t.isPrime(rounds: 10) {
primes.second[count2] = t
count2 += 1
}
}
 
return primes
}
 
 
let primes = pierpoint(n: 250)
 
print("First 50 Pierpoint primes of the first kind: \(Array(primes.first.prefix(50)))\n")
print("First 50 Pierpoint primes of the second kind: \(Array(primes.second.prefix(50)))")
print()
print("250th Pierpoint prime of the first kind: \(primes.first[249])")
print("250th Pierpoint prime of the second kind: \(primes.second[249])")</syntaxhighlight>
 
{{out}}
 
<pre>First 50 Pierpoint primes of the first kind: [2, 3, 5, 7, 13, 17, 19, 37, 73, 97, 109, 163, 193, 257, 433, 487, 577, 769, 1153, 1297, 1459, 2593, 2917, 3457, 3889, 10369, 12289, 17497, 18433, 39367, 52489, 65537, 139969, 147457, 209953, 331777, 472393, 629857, 746497, 786433, 839809, 995329, 1179649, 1492993, 1769473, 1990657, 2654209, 5038849, 5308417, 8503057]
 
First 50 Pierpoint primes of the second kind: [2, 3, 5, 7, 11, 17, 23, 31, 47, 53, 71, 107, 127, 191, 383, 431, 647, 863, 971, 1151, 2591, 4373, 6143, 6911, 8191, 8747, 13121, 15551, 23327, 27647, 62207, 73727, 131071, 139967, 165887, 294911, 314927, 442367, 472391, 497663, 524287, 786431, 995327, 1062881, 2519423, 10616831, 17915903, 18874367, 25509167, 30233087]
 
250th Pierpoint prime of the first kind: 62518864539857068333550694039553
250th Pierpoint prime of the second kind: 4111131172000956525894875083702271</pre>
 
=={{header|Wren}}==
{{trans|Go}}
{{libheader|Wren-big}}
{{libheader|Wren-fmt}}
The 3-smooth version. Just the first 250 - a tolerable 14 seconds or so on my machine.
<syntaxhighlight lang="wren">import "./big" for BigInt
import "./fmt" for Fmt
 
var pierpont = Fn.new { |n, first|
var p = [ List.filled(n, null), List.filled(n, null) ]
for (i in 0...n) {
p[0][i] = BigInt.zero
p[1][i] = BigInt.zero
}
p[0][0] = BigInt.two
var count = 0
var count1 = 1
var count2 = 0
var s = [BigInt.one]
var i2 = 0
var i3 = 0
var k = 1
while (count < n) {
var n2 = s[i2] * 2
var n3 = s[i3] * 3
var t
if (n2 < n3) {
t = n2
i2 = i2 + 1
} else {
t = n3
i3 = i3 + 1
}
if (t > s[k-1]) {
s.add(t.copy())
k = k + 1
t = t.inc
if (count1 < n && t.isProbablePrime(5)) {
p[0][count1] = t.copy()
count1 = count1 + 1
}
t = t - 2
if (count2 < n && t.isProbablePrime(5)) {
p[1][count2] = t.copy()
count2 = count2 + 1
}
count = count1.min(count2)
}
}
return p
}
 
var p = pierpont.call(250, true)
System.print("First 50 Pierpont primes of the first kind:")
for (i in 0...50) {
Fmt.write("$8i ", p[0][i])
if ((i-9)%10 == 0) System.print()
}
System.print("\nFirst 50 Pierpont primes of the second kind:")
for (i in 0...50) {
Fmt.write("$8i ", p[1][i])
if ((i-9)%10 == 0) System.print()
}
 
System.print("\n250th Pierpont prime of the first kind: %(p[0][249])")
System.print("\n250th Pierpont prime of the second kind: %(p[1][249])")</syntaxhighlight>
 
{{out}}
<pre>
First 50 Pierpont primes of the first kind:
2 3 5 7 13 17 19 37 73 97
109 163 193 257 433 487 577 769 1153 1297
1459 2593 2917 3457 3889 10369 12289 17497 18433 39367
52489 65537 139969 147457 209953 331777 472393 629857 746497 786433
839809 995329 1179649 1492993 1769473 1990657 2654209 5038849 5308417 8503057
 
First 50 Pierpont primes of the second kind:
2 3 5 7 11 17 23 31 47 53
71 107 127 191 383 431 647 863 971 1151
2591 4373 6143 6911 8191 8747 13121 15551 23327 27647
62207 73727 131071 139967 165887 294911 314927 442367 472391 497663
524287 786431 995327 1062881 2519423 10616831 17915903 18874367 25509167 30233087
 
250th Pierpont prime of the first kind: 62518864539857068333550694039553
 
250th Pierpont prime of the second kind: 4111131172000956525894875083702271
</pre>
 
=={{header|XPL0}}==
<syntaxhighlight lang "XPL0">func IsPrime(N); \Return 'true' if N is prime
int N, D;
[if N < 2 then return false;
if (N&1) = 0 then return N = 2;
if rem(N/3) = 0 then return N = 3;
D:= 5;
while D*D <= N do
[if rem(N/D) = 0 then return false;
D:= D+2;
if rem(N/D) = 0 then return false;
D:= D+4;
];
return true;
];
 
func Pierpont(N); \Return 'true' if N is a multiple of 2^U*3^V
int N;
[while (N&1) = 0 do N:= N>>1;
while N>1 do
[N:= N/3;
if rem(0) # 0 then return false;
];
return true;
];
 
int Kind, N, Count;
[Format(9, 0);
Kind:= 1;
repeat Count:= 0; N:= 2;
loop [if IsPrime(N) then
[if Pierpont(N-Kind) then
[Count:= Count+1;
RlOut(0, float(N));
if rem(Count/10) = 0 then CrLf(0);
if Count >= 50 then quit;
];
];
N:= N+1;
];
CrLf(0);
Kind:= -Kind;
until Kind = 1;
]</syntaxhighlight>
{{out}}
<pre>
2 3 5 7 13 17 19 37 73 97
109 163 193 257 433 487 577 769 1153 1297
1459 2593 2917 3457 3889 10369 12289 17497 18433 39367
52489 65537 139969 147457 209953 331777 472393 629857 746497 786433
839809 995329 1179649 1492993 1769473 1990657 2654209 5038849 5308417 8503057
 
2 3 5 7 11 17 23 31 47 53
71 107 127 191 383 431 647 863 971 1151
2591 4373 6143 6911 8191 8747 13121 15551 23327 27647
62207 73727 131071 139967 165887 294911 314927 442367 472391 497663
524287 786431 995327 1062881 2519423 10616831 17915903 18874367 25509167 30233087
 
</pre>
 
Line 1,301 ⟶ 3,081:
{{libheader|GMP}} GNU Multiple Precision Arithmetic Library
Using GMP's probabilistic primes makes it is easy and fast to test for primeness.
<langsyntaxhighlight lang="zkl">var [const] BI=Import("zklBigNum"); // libGMP
var [const] one=BI(1), two=BI(2), three=BI(3);
 
Line 1,329 ⟶ 3,109:
}
return(pps1,pps2)
}</langsyntaxhighlight>
<langsyntaxhighlight lang="zkl">pps1,pps2 := pierPonts(2_000);
 
println("The first 50 Pierpont primes (first kind):");
Line 1,341 ⟶ 3,121:
println("\n%4dth Pierpont prime, first kind: ".fmt(n), pps1[n-1]);
println( " second kind: ", pps2[n-1]);
}</langsyntaxhighlight>
{{out}}
<pre style="font-size:83%">
PureFox
9,476

edits

Retrieved from "https://rosettacode.org/wiki/Special:MobileDiff/293000...353893"