Pell numbers: Difference between revisions
→{{header|Oberon-07}}: Fix the near Isosoeles triangles
(Added Go) |
(→{{header|Oberon-07}}: Fix the near Isosoeles triangles) |
||
| (27 intermediate revisions by 14 users not shown) | |||
Line 94:
=={{header|ALGOL 68}}==
{{Trans|FreeBASIC|using Miller Rabin for finding Pell primes and some minor output format differences and showing 10 Pell primes}}
{{works with|ALGOL 68G|Any - tested with release 2.8.3.win32}}
Should work with any implementation of Algol 68 if LONG INT is large enough to hold the 90th Pell number/
{{libheader|ALGOL 68-primes}}
<syntaxhighlight lang="algol68">
BEGIN # find some Pell numbers - trans FreeBASIC ( which is trans Phix ) #
PR read "primes.incl.a68" PR
[ 0 : 90 ]LONG INT p, pl;
p[ 0 ] := 0; p[ 1 ] := 1;
pl[ 0 ] := 2; pl[ 1 ] := 2;
FOR n FROM 2 TO UPB p DO
p[ n ] := 2 * p[ n - 1 ] + p[ n - 2 ];
pl[ n ] := 2 * pl[ n - 1 ] + pl[ n - 2 ]
OD;
print( ( "First 20 Pell numbers:", newline ) );
FOR n FROM 0 TO 19 DO print( ( " ", whole( p[ n ], 0 ) ) ) OD;
print( ( newline, newline, "First 20 Pell-Lucas numbers:", newline ) );
FOR n FROM 0 TO 19 DO print( ( " ", whole( pl[ n ], 0 ) ) ) OD;
print( ( newline, newline, "First 20 rational approximations of sqrt(2) (" ) );
print( ( fixed( sqrt( 2 ), -15, 13 ), "):", newline ) );
FOR n TO 20 DO
LONG INT j = pl[ n ] OVER 2, d = p[ n ];
print( ( " ", whole( j, 0 ), "/", whole( d, 0 ), " ~= ", fixed( j / d, -15, 13 ), newline ) )
OD;
print( ( newline, "First 10 Pell primes:", newline, "index Pell prime", newline ) );
INT c := 0;
FOR pdx FROM 2 WHILE c < 10 DO
IF is probably prime( p[ pdx ] ) THEN
print( ( whole( pdx, -5 ), " ", whole( p[ pdx ], 0 ), newline ) );
c +:= 1
FI
OD;
print( ( newline, newline, "First 20 Newman-Shank-Williams numbers:", newline ) );
FOR n FROM 0 TO 19 DO
LONG INT nsw = p[ 2 * n ] + p[ 2 * n + 1 ];
print( ( " ", whole( nsw, 0 ) ) ); IF n = 13 THEN print( ( newline ) ) FI
OD;
print( ( newline, newline, "First 20 near isosceles right triangles:", newline ) );
LONG INT i0 := 0, i1 := 1, t := 1, found := 0;
FOR i FROM 2 WHILE found < 20 DO
LONG INT i2 = i1*2 + i0;
IF ODD i THEN
print( ( " [", whole( t, 0 ), ", ", whole( t + 1, 0 ), ", ", whole( i2, 0 ), "]", newline ) );
found +:= 1
FI;
t +:= i2; i0 := i1; i1 := i2
OD
END
</syntaxhighlight>
{{out}}
<pre style="height:40ex;overflow:scroll;">
First 20 Pell numbers:
0 1 2 5 12 29 70 169 408 985 2378 5741 13860 33461 80782 195025 470832 1136689 2744210 6625109
First 20 Pell-Lucas numbers:
2 2 6 14 34 82 198 478 1154 2786 6726 16238 39202 94642 228486 551614 1331714 3215042 7761798 18738638
First 20 rational approximations of sqrt(2) (1.4142135623731):
1/1 ~= 1.0000000000000
3/2 ~= 1.5000000000000
7/5 ~= 1.4000000000000
17/12 ~= 1.4166666666667
41/29 ~= 1.4137931034483
99/70 ~= 1.4142857142857
239/169 ~= 1.4142011834320
577/408 ~= 1.4142156862745
1393/985 ~= 1.4142131979695
3363/2378 ~= 1.4142136248949
8119/5741 ~= 1.4142135516461
19601/13860 ~= 1.4142135642136
47321/33461 ~= 1.4142135620573
114243/80782 ~= 1.4142135624273
275807/195025 ~= 1.4142135623638
665857/470832 ~= 1.4142135623747
1607521/1136689 ~= 1.4142135623728
3880899/2744210 ~= 1.4142135623731
9369319/6625109 ~= 1.4142135623731
22619537/15994428 ~= 1.4142135623731
First 10 Pell primes:
index Pell prime
2 2
3 5
5 29
11 5741
13 33461
29 44560482149
41 1746860020068409
53 68480406462161287469
59 13558774610046711780701
89 4125636888562548868221559797461449
First 20 Newman-Shank-Williams numbers:
1 7 41 239 1393 8119 47321 275807 1607521 9369319 54608393 318281039 1855077841 10812186007
63018038201 367296043199 2140758220993 12477253282759 72722761475561 423859315570607
First 20 near isosceles right triangles:
[3, 4, 5]
[20, 21, 29]
[119, 120, 169]
[696, 697, 985]
[4059, 4060, 5741]
[23660, 23661, 33461]
[137903, 137904, 195025]
[803760, 803761, 1136689]
[4684659, 4684660, 6625109]
[27304196, 27304197, 38613965]
[159140519, 159140520, 225058681]
[927538920, 927538921, 1311738121]
[5406093003, 5406093004, 7645370045]
[31509019100, 31509019101, 44560482149]
[183648021599, 183648021600, 259717522849]
[1070379110496, 1070379110497, 1513744654945]
[6238626641379, 6238626641380, 8822750406821]
[36361380737780, 36361380737781, 51422757785981]
[211929657785303, 211929657785304, 299713796309065]
[1235216565974040, 1235216565974041, 1746860020068409]
</pre>
=={{header|ALGOL W}}==
{{Trans|FreeBASIC|via Algol 68 - basic task only and only 5 Pell primes as Algol W limits integers to 32-bit.}}
<syntaxhighlight lang="algolw">
begin % find some Pell numbers - trans FreeBasic ( which is trans Phix ) %
% returns true if n is prime, false otherwise, uses trial division %
logical procedure isPrime ( integer value n ) ;
if n < 3 then n = 2
else if n rem 3 = 0 then n = 3
else if not odd( n ) then false
else begin
logical prime;
integer f, f2, toNext;
prime := true;
f := 5;
f2 := 25;
toNext := 24; % note: ( 2n + 1 )^2 - ( 2n - 1 )^2 = 8n %
while f2 <= n and prime do begin
prime := n rem f not = 0;
f := f + 2;
f2 := toNext;
toNext := toNext + 8
end while_f2_le_n_and_prime ;
prime
end isPrime ;
integer MAX_P;
MAX_P := 9;
begin
integer array p, pl ( 0 :: 20 ); % need more than 10 Pell numbers %
integer c, pdx; % to find the fifth Pell prime %
p( 0 ) := 0; p( 1 ) := 1;
pl( 0 ) := 2; pl( 1 ) := 2;
for n := 2 until 20 do begin
p( n ) := 2 * p( n - 1 ) + p( n - 2 );
pl( n ) := 2 * pl( n - 1 ) + pl( n - 2 )
end for_n ;
write( "First 10 Pell numbers:" );
for n := 0 until MAX_P do begin writeon( i_w := 1, s_w := 0, " ", p( n ) ) end;
write();write( "First 10 Pell-Lucas numbers:" );
for n := 0 until MAX_P do begin
writeon( i_w := 1, s_w := 0, " ", pl( n ) )
end for_n ;
write( s_w := 0, "First 10 rational approximations of sqrt(2) (" );
writeon( s_w := 0, r_format := "A", r_w := 8, r_d := 6, sqrt( 2 ), "):" );
for n := 1 until MAX_P do begin
integer j;
j := pl( n ) div 2;
write( i_w := 1, s_w := 0, " ", j, "/", p( n ), " ~= "
, r_format := "A", r_w := 8, r_d := 6, j / p( n )
)
end for_n;
write();write( "First 5 Pell primes:" ); write( "index Pell prime" );
c := 0;
pdx := 2;
while c < 5 do begin
if isPrime( p( pdx ) ) then begin
write( i_w := 5, s_w := 0, pdx, " ", i_w := 1, p( pdx ) );
c := c + 1
end if_isPrime_p_pdx ;
pdx := pdx + 1
end while_c_lt_5 ;
write(); write( "First 10 Newman-Shank-Williams numbers:" );write();
for n := 0 until MAX_P do begin
integer nsw;
nsw := p( 2 * n ) + p( 2 * n + 1 );
writeon( i_w := 1, s_w := 0, " ", nsw )
end for_n;
write();write( "First 10 near isosceles right triangles:" );
begin
integer i, i0, i1, i2, t, found;
i0 := 0; i1 := 1; t := 1; found := 0;
i := 1;
while found < 10 do begin
i := i + 1;
i2 := i1*2 + i0;
if odd( i ) then begin
write( i_w := 1, s_w := 0, " [", t, ", ", t + 1, ", ", i2, "]" );
found := found + 1
end if_odd_i ;
t := t + i2; i0 := i1; i1 := i2
end while_found_lt_10
end
end
end.
</syntaxhighlight>
{{out}}
<pre style="height:40ex;overflow:scroll;">
First 10 Pell numbers: 0 1 2 5 12 29 70 169 408 985
First 10 Pell-Lucas numbers: 2 2 6 14 34 82 198 478 1154 2786
First 10 rational approximations of sqrt(2) (1.414214):
1/1 ~= 1.000000
3/2 ~= 1.500000
7/5 ~= 1.400000
17/12 ~= 1.416667
41/29 ~= 1.413793
99/70 ~= 1.414286
239/169 ~= 1.414201
577/408 ~= 1.414216
1393/985 ~= 1.414213
First 5 Pell primes:
index Pell prime
2 2
3 5
5 29
11 5741
13 33461
First 10 Newman-Shank-Williams numbers:
1 7 41 239 1393 8119 47321 275807 1607521 9369319
First 10 near isosceles right triangles:
[3, 4, 5]
[20, 21, 29]
[119, 120, 169]
[696, 697, 985]
[4059, 4060, 5741]
[23660, 23661, 33461]
[137903, 137904, 195025]
[803760, 803761, 1136689]
[4684659, 4684660, 6625109]
[27304196, 27304197, 38613965]
</pre>
=={{header|Common Lisp}}==
<syntaxhighlight lang="CommonLisp">
(defun recurrent-sequence-2 (a0 a1 k1 k2 max)
"A generic function for any recurrent sequence of order 2, where a0 and a1 are the initial elements,
k1 is the factor of a(n-1) and k2 is the factor of a(n-2)"
(do* ((i 0 (1+ i))
(result (list a1 a0))
(b0 a0 b1)
(b1 a1 b2)
(b2 (+ (* k1 b1) (* k2 b0)) (+ (* k1 b1) (* k2 b0))) )
((> i max) (nreverse result))
(push b2 result) ))
(defun pell-sequence (max)
(recurrent-sequence-2 0 1 2 1 max) )
(defun pell-lucas-sequence (max)
(recurrent-sequence-2 2 2 2 1 max) )
(defun fibonacci-sequence (max) ; As an extra bonus, you get Fibonacci's numbers with this simple call
(recurrent-sequence-2 1 1 1 1 max) )
(defun rational-approximation-sqrt2 (max)
"Approximate square root of 2 with (P(n-1)+P(n-2))/P(n)"
(butlast (maplist #'(lambda (l) (/ (+ (first l) (or (second l) 0)) (or (second l) 1))) (pell-sequence max))) )
(defun pell-primes (max)
(do* ((i 0 (1+ i))
(result (list 1 0))
(indices nil)
(b0 0 b1)
(b1 1 b2)
(b2 (+ (* 2 b1) (* 1 b0)) (+ (* 2 b1) (* 1 b0))) )
((> (length result) max) (values (nreverse result)(nreverse indices)))
; primep can be any function determining whether a number is prime, for example,
; https://rosettacode.org/wiki/Primality_by_Wilson%27s_theorem#Common_Lisp
(when (primep b2)
(push b2 result)
(push i indices) )))
</syntaxhighlight>
The first 10 Pell numbers are:
<syntaxhighlight lang="CommonLisp">
(pell-sequence 10)
(0 1 2 5 12 29 70 169 408 985 2378 5741 13860)
</syntaxhighlight>
The first 10 Pell-Lucas numbers are:
<syntaxhighlight lang="CommonLisp">
(pell-lucas-sequence 10)
(2 2 6 14 34 82 198 478 1154 2786 6726 16238 39202)
</syntaxhighlight>
Rational approximation of square root of 2:
<syntaxhighlight lang="CommonLisp">
(rational-approximation-sqrt2 10)
(1 3/2 7/5 17/12 41/29 99/70 239/169 577/408 1393/985 3363/2378 8119/5741
19601/13860)
</syntaxhighlight>
The same in decimal form:
<syntaxhighlight lang="CommonLisp">
(mapcar #'float (rational-approximation-sqrt2 10))
(1.0 1.5 1.4 1.4166666 1.4137931 1.4142857 1.4142011 1.4142157 1.4142132
1.4142137 1.4142135 1.4142135)
</syntaxhighlight>
First 7 Pell primes and, as second value, the corresponding indices:
<syntaxhighlight lang="CommonLisp">
(pell-primes 7)
(0 1 2 5 29 5741 33461 44560482149) ;
(0 1 3 9 11 27)
</syntaxhighlight>
Still missing some of the taks
=={{header|FreeBASIC}}==
{{trans|Phix}}
<
Function isPrime(Byval ValorEval As Integer) As Boolean
Line 165 ⟶ 508:
i += 1
Loop
Sleep</
{{out}}
<pre style="height:40ex;overflow:scroll;">First 20 Pell numbers:
Line 235 ⟶ 578:
{{trans|Wren}}
{{libheader|Go-rcu}}
<
import (
Line 308 ⟶ 651:
u1 = u2
}
}</
{{out}}
Line 385 ⟶ 728:
(1235216565974040, 1235216565974041, 1746860020068409)
</pre>
=={{header|Haskell}}==
<syntaxhighlight lang="haskell">import Data.Numbers.Primes (isPrime)
----------------------- PELL SERIES ----------------------
pell :: Integer -> Integer -> [Integer]
pell a b = a : b : zipWith (+) (pell a b) ((2 *) <$> tail (pell a b))
a000129, a002203, a001333, a086383, a096650, a002315 :: [Integer]
a000129 = pell 0 1
a002203 = pell 2 2
a001333 = (`div` 2) <$> a002203
a086383 = filter isPrime a000129
a096650 = zip [0 ..] a000129 >>= (\(i, n) -> [i | isPrime n])
a002315 = 1 : 7 : zipWith (-) ((6 *) <$> tail a002315) a002315
------------------- PYTHAGOREAN TRIPLES ------------------
pythagoreanTriples :: [(Integer, Integer, Integer)]
pythagoreanTriples =
(tail . concat) $ zipWith3 go [0 ..] a000129 (scanl (+) 0 a000129)
where
go i p m
| odd i = [(m, succ m, p)]
| otherwise = []
-------------------------- TESTS -------------------------
main :: IO ()
main = do
mapM_
(\(k, xs) -> putStrLn ('\n' : k) >> print (take 10 xs))
[ ("a000129", a000129)
, ("a002203", a002203)
, ("a001333", a001333)
-- Waste of electrical power ?
-- ("a086383", a086383)
-- ("a096650", a096650
, ("a002315", a002315)
]
putStrLn "\nRational approximations to sqrt 2:"
mapM_ putStrLn $
(take 10 . tail) $
zipWith
(\n d ->
show n <>
('/' : show d) <> " -> " <> show (fromIntegral n / fromIntegral d))
a001333
a000129
putStrLn "\nPythagorean triples:"
mapM_ print $ take 10 pythagoreanTriples</syntaxhighlight>
{{Out}}
<pre>
a000129
[0,1,2,5,12,29,70,169,408,985]
a002203
[2,2,6,14,34,82,198,478,1154,2786]
a001333
[1,1,3,7,17,41,99,239,577,1393]
a002315
[1,7,41,239,1393,8119,47321,275807,1607521,9369319]
Rational approximations to sqrt 2:
1/1 -> 1.0
3/2 -> 1.5
7/5 -> 1.4
17/12 -> 1.4166666666666667
41/29 -> 1.4137931034482758
99/70 -> 1.4142857142857144
239/169 -> 1.4142011834319526
577/408 -> 1.4142156862745099
1393/985 -> 1.4142131979695431
3363/2378 -> 1.4142136248948696
Pythagorean triples:
(3,4,5)
(20,21,29)
(119,120,169)
(696,697,985)
(4059,4060,5741)
(23660,23661,33461)
(137903,137904,195025)
(803760,803761,1136689)
(4684659,4684660,6625109)
(27304196,27304197,38613965)</pre>
=={{header|J}}==
Line 391 ⟶ 830:
For example:
<
Pn=: (%:8) %~(1+%:2)&^ - (1-%:2)&^ NB. pell, closed form
Qn=: (1+%:2)&^ + (1-%:2)&^ NB. pell lucas, closed form
QN=: +: %&Pn ] NB. pell lucas, closed form
qn=: 2 * (+&Pn <:) NB. pell lucas, closed form</
Thus:
<
0 1 2 5 12 29 70 169 408 985 2378
Pn i.11
Line 410 ⟶ 849:
0 2 6 14 34 82 198 478 1154 2786 6726
qn i.11
2 2 6 14 34 82 198 478 1154 2786 6726</
QN (which is defined as P<sub>2n</sub>/P<sub>n</sub>) doesn't get the first element of the pell lucas sequence right. We could fix this by changing the definition:
<
QN i.11
2 2 6 14 34 82 198 478 1154 2786 6726</
Continuing... the first ten rational approximations to √2 here would be:
<
1 1.5 1.4 1.41667 1.41379 1.41429 1.4142 1.41422 1.41421 1.41421
}.(%~ _1}. +//.@,:~) nextPell^:9(0 1x)
1 3r2 7r5 17r12 41r29 99r70 239r169 577r408 1393r985 3363r2378</
The first ten pell primes are:
<
2 5 29 5741 33461 44560482149 1746860020068409 68480406462161287469 13558774610046711780701 4125636888562548868221559797461449</
Their indices are:
<
2 3 5 11 13 29 41 53 59 89</
The NSW numbers are the sums of (non-overlapping) pairs of pell numbers, or:
<
1 7 41 239 1393 8119 47321 275807 1607521 9369319 54608393
</syntaxhighlight>
The first ten pell based pythogorean triples would be:
<
3 4 5
20 21 29
Line 448 ⟶ 887:
803760 803761 1136689
4684659 4684660 6625109
27304196 27304197 38613965</
=={{header|Java}}==
<syntaxhighlight lang="java">
import java.math.BigDecimal;
import java.math.BigInteger;
import java.math.MathContext;
import java.util.ArrayList;
import java.util.List;
public final class PellNumbers {
public static void main(String[] aArgs) {
System.out.println("7 Tasks");
System.out.println("-------");
System.out.println("Task 1, Pell Numbers: " + pellNumbers(TERM_COUNT));
System.out.println();
System.out.println("Task 2, Pell-Lucas Numbers: " + pellLucasNumbers(TERM_COUNT));
System.out.println();
System.out.println("Task 3, Approximations to square root of 2:");
List<Rational> rationals = squareRoot2(TERM_COUNT + 1);
for ( int i = 1; i < TERM_COUNT + 1; i++ ) {
System.out.println(rationals.get(i) + " = " + rationals.get(i).toBigDecimal());
}
System.out.println();
List<Pair> pairs = pellPrimes(TERM_COUNT);
System.out.println("Task 4, Pell primes:");
for ( int i = 0; i < TERM_COUNT; i++ ) {
System.out.println(pairs.get(i).pellPrime);
}
System.out.println();
System.out.print("Task 5, Pell indices of Pell primes:");
for ( int i = 0; i < TERM_COUNT; i++ ) {
System.out.print(pairs.get(i).index + " ");
}
System.out.println();
System.out.println();
System.out.println("Task 6, Newman-Shank-Williams numbers: " + newmanShankWilliams(TERM_COUNT));
System.out.println();
System.out.println("Task , Near isoseles right triangles:");
List<Triple> nearIsoselesRightTriangles = nearIsoslesRightTriangles(TERM_COUNT);
for ( int i = 0; i < TERM_COUNT; i++ ) {
System.out.println(nearIsoselesRightTriangles.get(i));
}
}
private static List<BigInteger> pellNumbers(int aTermCount) {
PellIterator pellIterator = new PellIterator(BigInteger.ZERO, BigInteger.ONE);
List<BigInteger> result = new ArrayList<BigInteger>();
for ( int i = 0; i < aTermCount; i++ ) {
result.add(pellIterator.next());
}
return result;
}
private static List<BigInteger> pellLucasNumbers(int aTermCount) {
PellIterator pellLucasIterator = new PellIterator(BigInteger.TWO, BigInteger.TWO);
List<BigInteger> result = new ArrayList<BigInteger>();
for ( int i = 0; i < aTermCount; i++ ) {
result.add(pellLucasIterator.next());
}
return result;
}
private static List<Rational> squareRoot2(int aTermCount) {
PellIterator pellIterator = new PellIterator(BigInteger.ZERO, BigInteger.ONE);
PellIterator pellLucasIterator = new PellIterator(BigInteger.TWO, BigInteger.TWO);
List<Rational> result = new ArrayList<Rational>();
for ( int i = 0; i < aTermCount; i++ ) {
result.add( new Rational(pellLucasIterator.next().divide(BigInteger.TWO), pellIterator.next()) );
}
return result;
}
private static List<Pair> pellPrimes(int aTermCount) {
PellIterator pellIterator = new PellIterator(BigInteger.ZERO, BigInteger.ONE);
int index = 0;
int count = 0;
List<Pair> result = new ArrayList<Pair>();
while ( count < aTermCount ) {
BigInteger pellNumber = pellIterator.next();
if ( pellNumber.isProbablePrime(16) ) {
result.add( new Pair(pellNumber, index) );
count += 1;
}
index += 1;
}
return result;
}
private static List<BigInteger> newmanShankWilliams(int aTermCount) {
PellIterator pellIterator = new PellIterator(BigInteger.ZERO, BigInteger.ONE);
List<BigInteger> result = new ArrayList<BigInteger>();
for ( int i = 0; i < aTermCount; i++ ) {
BigInteger pellNumber = pellIterator.next();
result.add(pellNumber.add(pellIterator.next()));
}
return result;
}
private static List<Triple> nearIsoslesRightTriangles(int aTermCount) {
PellIterator pellIterator = new PellIterator(BigInteger.ZERO, BigInteger.ONE);
pellIterator.next();
List<Triple> result = new ArrayList<Triple>();
BigInteger sum = pellIterator.next();
for ( int i = 0; i < aTermCount; i++ ) {
sum = sum.add(pellIterator.next());
BigInteger nextTerm = pellIterator.next();
result.add( new Triple(sum, sum.add(BigInteger.ONE), nextTerm) );
sum = sum.add(nextTerm);
}
return result;
}
private static class PellIterator {
public PellIterator(BigInteger aFirst, BigInteger aSecond) {
a = aFirst; b = aSecond;
}
public BigInteger next() {
aCopy = a;
bCopy = b;
b = b.add(b).add(a);
a = bCopy;
return aCopy;
}
private BigInteger a, aCopy, b, bCopy;
}
private static record Rational(BigInteger numerator, BigInteger denominator) {
public BigDecimal toBigDecimal() {
return new BigDecimal(numerator).divide( new BigDecimal(denominator), mathContext );
}
public String toString() {
return numerator + " / " + denominator;
}
private static MathContext mathContext = new MathContext(34);
}
private static record Pair(BigInteger pellPrime, int index) {}
private static record Triple(BigInteger shortSide, BigInteger longSide, BigInteger hypotenuse) {
public String toString() {
return "(" + shortSide + ", " + longSide + ", " + hypotenuse + ")" ;
}
}
private static final int TERM_COUNT = 10;
}
</syntaxhighlight>
{{ out }}
<pre>
7 Tasks
-------
Task 1, Pell Numbers: [0, 1, 2, 5, 12, 29, 70, 169, 408, 985]
Task 2, Pell-Lucas Numbers: [2, 2, 6, 14, 34, 82, 198, 478, 1154, 2786]
Task 3, Approximations to square root of 2:
1 / 1 = 1
3 / 2 = 1.5
7 / 5 = 1.4
17 / 12 = 1.416666666666666666666666666666667
41 / 29 = 1.413793103448275862068965517241379
99 / 70 = 1.414285714285714285714285714285714
239 / 169 = 1.414201183431952662721893491124260
577 / 408 = 1.414215686274509803921568627450980
1393 / 985 = 1.414213197969543147208121827411168
3363 / 2378 = 1.414213624894869638351555929352397
Task 4, Pell primes:
2
5
29
5741
33461
44560482149
1746860020068409
68480406462161287469
13558774610046711780701
4125636888562548868221559797461449
Task 5, Pell indices of Pell primes:2 3 5 11 13 29 41 53 59 89
Task 6, Newman-Shank-Williams numbers: [1, 7, 41, 239, 1393, 8119, 47321, 275807, 1607521, 9369319]
Task , Near isoseles right triangles:
(3, 4, 5)
(20, 21, 29)
(119, 120, 169)
(696, 697, 985)
(4059, 4060, 5741)
(23660, 23661, 33461)
(137903, 137904, 195025)
(803760, 803761, 1136689)
(4684659, 4684660, 6625109)
(27304196, 27304197, 38613965)
</pre>
=={{header|Julia}}==
<
function pellnumbers(wanted)
Line 517 ⟶ 1,163:
printrows("Twenty Newman-Shank-Williams numbers:", newmanshankwilliams(20), 17, 5)
printrows("Twenty near isosceles triangle triplets:", nearisosceles(20), 52, 2)
</
<pre>
Twenty Pell numbers:
Line 577 ⟶ 1,223:
(211929657785303, 211929657785304, 299713796309065) (1235216565974040, 1235216565974041, 1746860020068409)
</pre>
=={{header|Mathematica}}/{{header|Wolfram Language}}==
<syntaxhighlight lang="mathematica">ClearAll[PellNumber, PellLucasNumber]
PellNumber[0] = 0;
PellNumber[1] = 1;
PellNumber[n_] := PellNumber[n] = 2 PellNumber[n - 1] + PellNumber[n - 2]
PellLucasNumber[0] = 2;
PellLucasNumber[1] = 2;
PellLucasNumber[n_] := PellLucasNumber[n] = 2 PellLucasNumber[n - 1] + PellLucasNumber[n - 2]
pns = PellNumber /@ Range[0, 9]
plns = PellLucasNumber /@ Range[0, 9]
den = Rest@pns;
num = Rest@plns/2;
approx = num/den
N[approx]
pns = {#, PellNumber[#]} & /@ Range[0, 100];
Select[pns, Last/*PrimeQ, 10] // Grid
ClearAll[PellS]
PellS[n_] := If[n == 0, 1, PellNumber[2 n] + PellNumber[2 n + 1]]
PellS /@ Range[0, 19]
ClearAll[PythagoreanTriple]
PythagoreanTriple[n_Integer] := Module[{hypo, short, long},
hypo = PellNumber[2 n + 1];
short = Total[PellNumber /@ Range[2 n]];
long = short + 1;
{short, long, hypo}
]
PythagoreanTriple /@ Range[10]</syntaxhighlight>
{{out}}
<pre>{0, 1, 2, 5, 12, 29, 70, 169, 408, 985}
{2, 2, 6, 14, 34, 82, 198, 478, 1154, 2786}
{1, 3/2, 7/5, 17/12, 41/29, 99/70, 239/169, 577/408, 1393/985}
{1., 1.5, 1.4, 1.41667, 1.41379, 1.41429, 1.4142, 1.41422, 1.41421}
2 2
3 5
5 29
11 5741
13 33461
29 44560482149
41 1746860020068409
53 68480406462161287469
59 13558774610046711780701
89 4125636888562548868221559797461449
{1, 7, 41, 239, 1393, 8119, 47321, 275807, 1607521, 9369319, 54608393, 318281039, 1855077841, 10812186007, 63018038201, 367296043199, 2140758220993, 12477253282759, 72722761475561, 423859315570607}
{{3, 4, 5}, {20, 21, 29}, {119, 120, 169}, {696, 697, 985}, {4059, 4060, 5741}, {23660, 23661, 33461}, {137903, 137904, 195025}, {803760, 803761, 1136689}, {4684659, 4684660, 6625109}, {27304196, 27304197, 38613965}}</pre>
=={{header|Nim}}==
<syntaxhighlight lang="Nim">import std/[strformat, strutils, sugar]
import integers
template isOdd(n: int): bool = (n and 1) != 0
iterator pellSequence[T](first, second: T; lim = -1): (int, T) =
## Yield the sucessive values of a Pell or Pell-Lucas sequence
## preceded by their rank.
## If "lim" is specified and greater than 2, only the "lim" first
## values are computed.
## The iterator works with "int" values or "Integer" values.
var prev = first
var curr = second
yield (0, prev)
yield (1, curr)
var count = 2
while true:
swap prev, curr
curr += 2 * prev
yield (count, curr)
inc count
if count == lim: break
echo "First 10 Pell numbers:"
let p = collect:
for (idx, val) in pellSequence(0, 1, 11): val
echo p[0..9].join(" ")
echo "\nFirst 10 Pell-Lucas numbers:"
let q = collect:
for (idx, val) in pellSequence(2, 2, 11): val
echo q[0..9].join(" ")
echo "\nFirst 10 rational approximations of √2:"
for i in 1..10:
let n = q[i] div 2
let d = p[i]
let r = &"{n}/{d}"
echo &"{r:>9} = {n/d:.17}"
echo "\nFirst 10 Pell primes:"
# To avoid an overflow, we need to use Integer values here.
var indices: seq[int]
var count = 0
for (idx, p) in pellSequence(newInteger(0), newInteger(1)):
if p.isPrime:
echo p
indices.add idx
inc count
if count == 10: break
echo "\nFirst 10 Pell primes indices:"
echo indices.join(" ")
echo "\nFirst 10 Newman-Shank-Williams numbers:"
count = 0
var prev: int
for (idx, p) in pellSequence(0, 1):
if idx.isOdd:
inc count
stdout.write prev + p
if count == 10: break
stdout.write ' '
else:
prev = p
echo()
echo "\nFirst 10 near isosceles right triangles:"
count = 0
var sum = 0
for (idx, p) in pellSequence(0, 1):
if idx.isOdd and sum != 0:
echo (sum, sum + 1, p)
inc count
if count == 10: break
inc sum, p
</syntaxhighlight>
{{out}}
<pre>First 10 Pell numbers:
0 1 2 5 12 29 70 169 408 985
First 10 Pell-Lucas numbers:
2 2 6 14 34 82 198 478 1154 2786
First 10 rational approximations of √2:
1/1 = 1.0000000000000000
3/2 = 1.5000000000000000
7/5 = 1.3999999999999999
17/12 = 1.4166666666666667
41/29 = 1.4137931034482758
99/70 = 1.4142857142857144
239/169 = 1.4142011834319526
577/408 = 1.4142156862745099
1393/985 = 1.4142131979695431
3363/2378 = 1.4142136248948696
First 10 Pell primes:
2
5
29
5741
33461
44560482149
1746860020068409
68480406462161287469
13558774610046711780701
4125636888562548868221559797461449
First 10 Pell primes indices:
2 3 5 11 13 29 41 53 59 89
First 10 Newman-Shank-Williams numbers:
1 7 41 239 1393 8119 47321 275807 1607521 9369319
First 10 near isosceles right triangles:
(3, 4, 5)
(20, 21, 29)
(119, 120, 169)
(696, 697, 985)
(4059, 4060, 5741)
(23660, 23661, 33461)
(137903, 137904, 195025)
(803760, 803761, 1136689)
(4684659, 4684660, 6625109)
(27304196, 27304197, 38613965)
</pre>
=={{header|Oberon-07}}==
{{Trans|ALGOL_W|which is a translation of FreeBASIC which is a translation of Phix}}
<syntaxhighlight lang="modula2">
MODULE PellNumbers; (* find some Pell numbersl trans Algol W, FreeBasic, Phix *)
IMPORT Out, Math;
CONST MAXP = 9;
VAR p, pl :ARRAY 21 OF INTEGER; (* need more than 10 Pell numbers *)
c, pdx, j, n, nsw :INTEGER; (* to find the fifth Pell prime *)
(* returns true if n is prime, false otherwise, uses trial division *)
PROCEDURE isPrime( n : INTEGER ):BOOLEAN;
VAR prime :BOOLEAN;
f, f2, toNext :INTEGER;
BEGIN
IF n < 3 THEN prime := n = 2
ELSIF n MOD 3 = 0 THEN prime := n = 3
ELSIF ~ ODD( n ) THEN prime := FALSE
ELSE
prime := TRUE;
f := 5;
f2 := 25;
toNext := 24; (* note: ( 2n + 1 )^2 - ( 2n - 1 )^2 = 8n *)
WHILE ( f2 <= n ) & prime DO
prime := n MOD f # 0;
f := f + 2;
f2 := toNext;
toNext := toNext + 8
END
END
RETURN prime
END isPrime;
PROCEDURE NearIsoscelesRightTriangles;
VAR i, i0, i1, i2, t, found :INTEGER;
BEGIN
i0 := 0; i1 := 1; t := 1; found := 0;
i := 1;
WHILE found < 10 DO
i := i + 1;
i2 := i1*2 + i0;
IF ODD( i ) THEN
Out.String( " [" );Out.Int( t, 0 );
Out.String( ", " );Out.Int( t + 1, 0 );
Out.String( ", " );Out.Int( i2, 0 );
Out.String( "]" );Out.Ln;
found := found + 1
END;
t := t + i2; i0 := i1; i1 := i2
END
END NearIsoscelesRightTriangles;
BEGIN
p[ 0 ] := 0; p[ 1 ] := 1;
pl[ 0 ] := 2; pl[ 1 ] := 2;
FOR n := 2 TO 20 DO
p[ n ] := 2 * p[ n - 1 ] + p[ n - 2 ];
pl[ n ] := 2 * pl[ n - 1 ] + pl[ n - 2 ]
END;
Out.String( "First 10 Pell numbers:" );
FOR n := 0 TO MAXP DO Out.String( " " );Out.Int( p[ n ], 1 ) END; Out.Ln;
Out.String( "First 10 Pell-Lucas numbers:" );
FOR n := 0 TO MAXP DO Out.String( " " );Out.Int( pl[ n ], 1 ) END; Out.Ln;
Out.Ln;Out.String( "First 10 rational approximations of sqrt(2) (" );
Out.Real( Math.sqrt( 2.0 ), 8 );Out.String( "):" );
FOR n := 1 TO MAXP DO
j := pl[ n ] DIV 2;
Out.Ln;
Out.String( " " );Out.Int( j, 0 );Out.String( "/" );Out.Int( p[ n ], 0 );
Out.String( " ~= " );Out.Real( FLT( j ) / FLT( p[ n ] ), 8 );
END;
Out.Ln;
Out.Ln;Out.String( "First 5 Pell primes:" );
Out.Ln;Out.String( "index Pell prime" );
c := 0;
pdx := 2;
WHILE c < 5 DO
IF isPrime( p[ pdx ] ) THEN
Out.Ln;Out.Int( pdx, 6 );Out.String( " " );Out.Int( p[ pdx ], 0 );
c := c + 1
END;
pdx := pdx + 1
END;
Out.Ln;Out.Ln;Out.String( "First 10 Newman-Shank-Williams numbers:" );Out.Ln;
FOR n := 0 TO MAXP DO
nsw := p[ 2 * n ] + p[ 2 * n + 1 ];
Out.String( " " );Out.Int( nsw, 0 )
END;
Out.Ln;
Out.Ln;Out.String( "First 10 near isosceles right triangles:" );Out.Ln;
NearIsoscelesRightTriangles
END PellNumbers.
</syntaxhighlight>
{{out}}
<pre style="height:40ex;overflow:scroll;">
First 10 Pell numbers: 0 1 2 5 12 29 70 169 408 985
First 10 Pell-Lucas numbers: 2 2 6 14 34 82 198 478 1154 2786
First 10 rational approximations of sqrt(2) (1.414214):
1/1 ~= 1.000000
3/2 ~= 1.500000
7/5 ~= 1.400000
17/12 ~= 1.416667
41/29 ~= 1.413793
99/70 ~= 1.414286
239/169 ~= 1.414201
577/408 ~= 1.414216
1393/985 ~= 1.414213
First 5 Pell primes:
index Pell prime
2 2
3 5
5 29
11 5741
13 33461
First 10 Newman-Shank-Williams numbers:
1 7 41 239 1393 8119 47321 275807 1607521 9369319
First 10 near isosceles right triangles:
[3, 4, 5]
[20, 21, 29]
[119, 120, 169]
[696, 697, 985]
[4059, 4060, 5741]
[23660, 23661, 33461]
[137903, 137904, 195025]
[803760, 803761, 1136689]
[4684659, 4684660, 6625109]
[27304196, 27304197, 38613965]
</pre>
=={{header|Perl}}==
{{libheader|ntheory}}
<syntaxhighlight lang="perl">use strict;
use warnings;
use feature <state say>;
use bignum;
use ntheory 'is_prime';
use List::Util <sum head>;
use List::Lazy 'lazy_list';
my $upto = 20;
my $pell_recur = lazy_list { state @p = (0, 1); push @p, 2*$p[1] + $p[0]; shift @p };
my $pell_lucas_recur = lazy_list { state @p = (2, 2); push @p, 2*$p[1] + $p[0]; shift @p };
my @pell;
push @pell, $pell_recur->next() for 1..1500; #wart;
my @pell_lucas;
push @pell_lucas, $pell_lucas_recur->next() for 1..$upto+1;
say "First $upto Pell numbers:";
say join ' ', @pell[0..$upto-1];
say "\nFirst $upto Pell-Lucas numbers:";
say join ' ', @pell_lucas[0..$upto-1];
say "\nFirst $upto rational approximations of √2:";
say sprintf "%d/%d - %1.16f", $pell[$_-1] + $pell[$_], $pell[$_], ($pell[$_-1]+$pell[$_])/$pell[$_] for 1..$upto;
say "\nFirst $upto Pell primes:";
say join "\n", head $upto, grep { is_prime $_ } @pell;
say "\nIndices of first $upto Pell primes:";
say join ' ', head $upto, grep { is_prime($pell[$_]) and $_ } 0..$#pell;
say "\nFirst $upto Newman-Shank-Williams numbers:";
say join ' ', map { $pell[2 * $_] + $pell[2 * $_+1] } 0..$upto-1;
say "\nFirst $upto near isoceles right tringles:";
map {
my $y = 2*$_ + 1;
my $x = sum @pell[0..$y-1];
printf "(%d, %d, %d)\n", $x, $x+1, $pell[$y]
} 1..$upto;</syntaxhighlight>
{{out}}
<pre style="height:20ex">First 20 Pell numbers:
0 1 2 5 12 29 70 169 408 985 2378 5741 13860 33461 80782 195025 470832 1136689 2744210 6625109
First 20 Pell-Lucas numbers:
2 2 6 14 34 82 198 478 1154 2786 6726 16238 39202 94642 228486 551614 1331714 3215042 7761798 18738638
First 20 rational approximations of √2:
1/1 - 1.0000000000000000
3/2 - 1.5000000000000000
7/5 - 1.4000000000000000
17/12 - 1.4166666666666667
41/29 - 1.4137931034482758
99/70 - 1.4142857142857144
239/169 - 1.4142011834319526
577/408 - 1.4142156862745099
1393/985 - 1.4142131979695431
3363/2378 - 1.4142136248948696
8119/5741 - 1.4142135516460548
19601/13860 - 1.4142135642135643
47321/33461 - 1.4142135620573204
114243/80782 - 1.4142135624272734
275807/195025 - 1.4142135623637995
665857/470832 - 1.4142135623746899
1607521/1136689 - 1.4142135623728214
3880899/2744210 - 1.4142135623731420
9369319/6625109 - 1.4142135623730870
22619537/15994428 - 1.4142135623730965
First 20 Pell primes:
2
5
29
5741
33461
44560482149
1746860020068409
68480406462161287469
13558774610046711780701
4125636888562548868221559797461449
4760981394323203445293052612223893281
161733217200188571081311986634082331709
2964793555272799671946653940160950323792169332712780937764687561
677413820257085084326543915514677342490435733542987756429585398537901
4556285254333448771505063529048046595645004014152457191808671945330235841
54971607658948646301386783144964782698772613513307493180078896702918825851648683235325858118170150873214978343601463118106546653220435805362395962991295556488036606954237309847762149971207793263738989
14030291214037674827921599320400561033992948898216351802670122530401263880575255235196727095109669287799074570417579539629351231775861429098849146880746524269269235328805333087546933690012894630670427794266440579064751300508834822795162874147983974059159392260220762973563561382652223360667198516093199367134903695783143116067743023134509886357032327271649
2434804314652199381956027075145741187716221548707931096877274520825143228915116227412484991366386864484767844200542482630246332092069382947111767723898168035847078557798454111405556629400142434835890123610082763986456199467423944182141028870863302603437534363208996458153115358483747994095302552907353919742211197822912892578751357668345638404394626711701120567186348490247426710813709165801137112237291901437566040249805155494297005186344325519103590369653438042689
346434895614929444828445967916634653215454504812454865104089892164276080684080254746939261017687341632569935171059945916359539268094914543114024020158787741692287531903178502306292484033576487391159597130834863729261484555671037916432206867189514675750227327687799973497042239286045783392065227614939379139866240959756584073664244580698830046194724340448293320938108876004367449471918175071251610962540447986139876845105399212429593945098472125140242905536711601925585608153109062121115635939797709
32074710952523740376423283403256578238321646122759160107427497117576305397686814013623874765833543023397971470911301264845142006214276865917420065183527313421909784286074786922242104480428021290764613639424408361555091057197776876849282654018358993099016644054242247557103410808928387071991436781136646322261169941417916607548507224950058710729258466238995253184617782314756913932650536663800753256087990078866003788647079369825102832504351225446531057648755795494571534144773842019836572551455718577614678081652481281009
Indices of first 20 Pell primes:
2 3 5 11 13 29 41 53 59 89 97 101 167 181 191 523 929 1217 1301 1361
First 20 Newman-Shank-Williams numbers:
1 7 41 239 1393 8119 47321 275807 1607521 9369319 54608393 318281039 1855077841 10812186007 63018038201 367296043199 2140758220993 12477253282759 72722761475561 423859315570607
First 20 near isosceles right triangles:
(3, 4, 5)
(20, 21, 29)
(119, 120, 169)
(696, 697, 985)
(4059, 4060, 5741)
(23660, 23661, 33461)
(137903, 137904, 195025)
(803760, 803761, 1136689)
(4684659, 4684660, 6625109)
(27304196, 27304197, 38613965)
(159140519, 159140520, 225058681)
(927538920, 927538921, 1311738121)
(5406093003, 5406093004, 7645370045)
(31509019100, 31509019101, 44560482149)
(183648021599, 183648021600, 259717522849)
(1070379110496, 1070379110497, 1513744654945)
(6238626641379, 6238626641380, 8822750406821)
(36361380737780, 36361380737781, 51422757785981)
(211929657785303, 211929657785304, 299713796309065)
(1235216565974040, 1235216565974041, 1746860020068409)</pre>
=={{header|Phix}}==
<!--<
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span>
<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: #000000;">0</span><span style="color: #0000FF;">,</span><span style="color: #000000;">1</span><span style="color: #0000FF;">},</span>
Line 620 ⟶ 1,717:
<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;">"[%d, %d, %d]\n"</span><span style="color: #0000FF;">,</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">side</span><span style="color: #0000FF;">,</span><span style="color: #000000;">side</span><span style="color: #0000FF;">+</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #000000;">hypot</span><span style="color: #0000FF;">})</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<!--</
{{out}}
<pre>
Line 696 ⟶ 1,793:
[211929657785303, 211929657785304, 299713796309065]
[1235216565974040, 1235216565974041, 1746860020068409]
</pre>
=={{header|Python}}==
<syntaxhighlight lang="python3">
# pell_numbers.py by Xing216
def is_prime(n):
if n == 1:
return False
i = 2
while i*i <= n:
if n % i == 0:
return False
i += 1
return True
def pell(p0: int,p1: int,its: int):
nums = [p0,p1]
primes = {}
idx = 2
while len(nums) != its:
p = 2*nums[-1]+nums[-2]
if is_prime(p):
primes[idx] = p
nums.append(p)
idx += 1
return nums, primes
def nsw(its: int,pell_nos: list):
nums = []
for i in range(its):
nums.append(pell_nos[2*i] + pell_nos[2*i+1])
return nums
def pt(its: int, pell_nos: list):
nums = []
for i in range(1,its+1):
hypot = pell_nos[2*i+1]
shorter_leg = sum(pell_nos[:2*i+1])
longer_leg = shorter_leg + 1
nums.append((shorter_leg,longer_leg,hypot))
return nums
pell_nos, pell_primes = pell(0,1,50)
pell_lucas_nos, pl_primes = pell(2,2,10)
nsw_nos = nsw(10, pell_nos)
pythag_triples = pt(10, pell_nos)
sqrt2_approx = {}
for idx, pell_no in enumerate(pell_nos):
numer = pell_nos[idx-1] + pell_no
if pell_no != 0:
sqrt2_approx[f"{numer}/{pell_no}"] = numer/pell_no
print(f"The first 10 Pell Numbers:\n {' '.join([str(_) for _ in pell_nos[:10]])}")
print(f"The first 10 Pell-Lucas Numbers:\n {' '.join([str(_) for _ in pell_lucas_nos])}")
print(f"The first 10 rational and decimal approximations of sqrt(2) ({(2**0.5):.10f}):")
print(" rational | decimal")
for rational in list(sqrt2_approx.keys())[:10]:
print(f"{rational:>10} ≈ {sqrt2_approx[rational]:.10f}")
print("The first 7 Pell Primes:")
print(" index | Pell Prime")
for idx, prime in pell_primes.items():
print(f"{idx:>6} | {prime}")
print(f"The first 10 Newman-Shank-Williams numbers:\n {' '.join([str(_) for _ in nsw_nos])}")
print(f"The first 10 near isosceles right triangles:")
for i in pythag_triples:
print(f" {i}")
</syntaxhighlight>
{{out}}
<pre style="height: 15em;">
The first 10 Pell Numbers:
0 1 2 5 12 29 70 169 408 985
The first 10 Pell-Lucas Numbers:
2 2 6 14 34 82 198 478 1154 2786
The first 10 rational and decimal approximations of sqrt(2) (1.4142135624):
rational | decimal
1/1 ≈ 1.0000000000
3/2 ≈ 1.5000000000
7/5 ≈ 1.4000000000
17/12 ≈ 1.4166666667
41/29 ≈ 1.4137931034
99/70 ≈ 1.4142857143
239/169 ≈ 1.4142011834
577/408 ≈ 1.4142156863
1393/985 ≈ 1.4142131980
3363/2378 ≈ 1.4142136249
The first 7 Pell Primes:
index | Pell Prime
2 | 2
3 | 5
5 | 29
11 | 5741
13 | 33461
29 | 44560482149
41 | 1746860020068409
The first 10 Newman-Shank-Williams numbers:
1 7 41 239 1393 8119 47321 275807 1607521 9369319
The first 10 near isosceles right triangles:
(3, 4, 5)
(20, 21, 29)
(119, 120, 169)
(696, 697, 985)
(4059, 4060, 5741)
(23660, 23661, 33461)
(137903, 137904, 195025)
(803760, 803761, 1136689)
(4684659, 4684660, 6625109)
(27304196, 27304197, 38613965)
</pre>
=={{header|Quackery}}==
<code>isprime</code> is defined at [[Primality by trial division#Quackery]].
As this method of testing for primality is not well suited to the task the Pell Primes part is limited to finding the first seven and their indices. More than that would be impractical.
<syntaxhighlight lang="Quackery"> [ $ "bigrat.qky" loadfile ] now!
[ ' [ 0 ] swap
witheach
[ over -1 peek
+ join ]
behead drop ] is cumsum ( [ --> [ )
[ dup -1 peek 2 *
over -2 peek + join ] is nextterm ( [ --> [ )
[ over 2 - times
nextterm
swap split drop ] is sequence ( n [ --> [ )
[ ' [ 0 1 ] sequence ] is pells ( n --> [ )
[ ' [ 2 2 ] sequence ] is companions ( n --> [ )
[ [] swap 1+
dup companions
behead drop
swap pells
behead drop
witheach
[ dip
[ behead 2 / ]
join nested
rot swap join
swap ]
drop ] is rootytwos ( n --> [ )
[ stack ] is index ( --> s )
[ temp put
1 index put
[] ' [ 0 1 ]
[ 1 index tally
over size
temp share < while
nextterm
behead drop
index share
isprime until
dup -1 peek
dup isprime iff
[ swap dip
[ index share
swap join
nested join ] ]
else drop
again ]
drop
index release
temp release ] is pellprimes ( n --> [ )
[ [] over 2 * pells
rot times
[ behead dip behead +
join ]
nip ] is nsws ( n --> [ )
[ [] swap 1+ dup
2 * 1+ pells
dup cumsum
swap rot times
[ over i^ 2 * peek
dup 1+ join
over i^ 2 * 1+ peek
join
dip rot nested join
unrot ]
2drop behead drop ] is nirts ( n --> [ )
say "Pell numbers " 10 pells echo
cr cr
say "Pell-Lucas's " 10 companions echo
cr cr
say "Approximations of sqrt(2)"
cr
10 rootytwos
witheach
[ do 2dup
vulgar$ echo$ sp
10 point$ echo$ cr ]
cr
say "Pell Primes "
7 pellprimes dup
[] swap
witheach [ 1 peek join ] echo
cr
say "their indices "
[] swap
witheach [ 0 peek join ] echo
cr cr
say "NSW numbers " 10 nsws echo
cr cr
say "Near isosceles right triangles"
cr
10 nirts witheach [ echo cr ]</syntaxhighlight>
{{out}}
<pre>Pell numbers [ 0 1 2 5 12 29 70 169 408 985 ]
Pell-Lucas's [ 2 2 6 14 34 82 198 478 1154 2786 ]
Approximations of sqrt(2)
1/1 1
3/2 1.5
7/5 1.4
17/12 1.4166666667
41/29 1.4137931034
99/70 1.4142857143
239/169 1.4142011834
577/408 1.4142156863
1393/985 1.414213198
3363/2378 1.4142136249
Pell Primes [ 2 5 29 5741 33461 44560482149 1746860020068409 ]
their indices [ 2 3 5 11 13 29 41 ]
NSW numbers [ 1 7 41 239 1393 8119 47321 275807 1607521 9369319 ]
Near isosceles right triangles
[ 3 4 5 ]
[ 20 21 29 ]
[ 119 120 169 ]
[ 696 697 985 ]
[ 4059 4060 5741 ]
[ 23660 23661 33461 ]
[ 137903 137904 195025 ]
[ 803760 803761 1136689 ]
[ 4684659 4684660 6625109 ]
[ 27304196 27304197 38613965 ]
</pre>
=={{header|Raku}}==
<syntaxhighlight lang="raku"
my $pell-lucas = lazy 2, 2, * + * × 2 … *;
Line 718 ⟶ 2,059:
say "\nFirst $upto near isosceles right triangles:";
map -> \p { printf "(%d, %d, %d)\n", |($_, $_+1 given $pell[^(2 × p + 1)].sum), $pell[2 × p + 1] }, 1..$upto;</
{{out}}
<pre style="height:40ex;overflow:scroll;">First 20 Pell numbers:
Line 797 ⟶ 2,138:
(211929657785303, 211929657785304, 299713796309065)
(1235216565974040, 1235216565974041, 1746860020068409)</pre>
=={{header|RPL}}==
Let's start first with a universal Pell sequence generator, which takes a smaller sequence as input:
{{works with|Halcyon Calc|4.2.8}}
{| class="wikitable"
! RPL code
! Comment
|-
|
≪ → n
≪ '''IF''' DUP SIZE n < '''THEN'''
LIST→ n SWAP -
1 SWAP '''START''' DUP2 DUP + + '''NEXT'''
n →LIST '''END'''
≫ ≫ ''''PELLS'''' STO
|
'''PELLS''' ''( { P(0)..P(m) } n -- { P(0)..P(n-1) } )''
if m < n
put P(0)..P(m) in the stack
loop n-m times: P(j+1) = 2*P(j) + P(j-1)
pop stack to a list
return new list
|}
We can then generate Pell and Pell-Lucas numbers with the same function:
{ 0 1 } 25 '''PELLS'''
{ 2 2 } 25 '''PELLS'''
{{out}}
<pre>
2: { 0 1 2 5 12 29 70 169 408 985 2378 5741 13860 33461 80782 195025 470832 1136689 2744210 6625109 15994428 38613965 93222358 225058681 543339720 }
1: { 2 2 6 14 34 82 198 478 1154 2786 6726 16238 39202 94642 228486 551614 1331714 3215042 7761798 18738638 45239074 109216786 263672646 636562078 1536796802 }
</pre>
Generating a sequence of rational approximations needs some code for appropriate formatting:
≪ { 0 1 } 16 '''PELLS''' → p
≪ { } 2 15 '''FOR''' j
p j GET "'" OVER p j 1 - GET + →STR + "/" + SWAP →STR + STR→ DUP →NUM ROT ROT + SWAP + '''NEXT'''
≫ ≫ ''''TASK3'''' STO
{{out}}
<pre>
{ '1/1' 1 '3/2' 1.5 '7/5' 1.4 '17/12' 1.41666666667 '41/29' 1.41379310345 '99/70' 1.41428571429 '239/169' 1.41420118343 '577/408' 1.41421568627 '1393/985' 1.41421319797 '3363/2378' 1.41421362489 '8119/5741' 1.41421355165 '19601/13860' 1.41421356421 '47321/33461' 1.41421356206 '114243/80782' 1.41421356243 }
</pre>
To look for prime Pell numbers, we need a fast Pell number generator and also a prime number checker, available [[Primality by trial division#RPL|here]]:
≪ '''DO'''
1 + ROT ROT DUP DUP + ROT + ROT
'''UNTIL''' 3 PICK '''BPRIM? END'''
≫ ≫ ''''PRPL'''' STO
≪ n { } → n primpell
≪ # 0d # 1d 1 1 n '''START '''
'''PRPL''' 3 PICK B→R OVER 1 - R→C primpell SWAP + 'primpell' STO '''NEXT'''
3 DROPN primpell
≫ ≫ ''''TASK4'''' STO
<code>'''TASK4'''</code> returns a list of complex numbers, the real part being the prime number and the imaginary part being the indice.
{{out}}
<pre>
1: { (2,2) (5,3) (29,5) (5741,11) (33461,13) (44560482149,29) }
</pre>
Unfortunately, the emulator's watchdog timer has refused to let us look beyond P<sub>29</sub>.
Generating the NSW sequence is again light work:
≪ → n
≪ { 0 1 } 2 n * 1 + '''PELLS''' { }
0 n 1 - '''FOR''' j
OVER j DUP + 1 + GET LAST 1 + GET + + '''NEXT'''
SWAP DROP
≫ ≫ ''''NSW'''' STO
15 '''NSW'''
{{out}}
<pre>
1: { 1 7 41 239 1393 8119 47321 275807 1607521 9369319 54608393 318281039 1855077841 10812186007 63018038201 }
</pre>
For dessert, the Pythagorean triplets can be found with:
≪ → n
≪ { 0 1 } 2 n * 2 + '''PELLS''' { }
1 n '''FOR''' j
OVER 2 j * 2 + GET LASTARG 1 - 1 SWAP SUB ∑LIST DUP 1 +
ROT 3 →LIST 1 →LIST + '''NEXT'''
SWAP DROP
≫ ≫ ''''PYTH3'''' STO
15 '''PYTH3'''
This program runs on HP-48G or newer models only, but it can be adapted for older ones by creating a homemade version of <code>∑LIST</code>, which sums the items of a list, and replacing <code>LASTARG</code> with <code>LAST</code>.
{{out}}
<pre>
1: { { 3 4 5 } { 20 21 29 } { 119 120 169 } { 696 697 985 } { 4059 4060 5741 } { 23660 23661 33461 } { 137903 137904 195025 } { 803760 803761 1136689 } { 4684659 4684660 6625109 } { 27304196 27304197 38613965 } { 159140519 159140520 225058681 } { 927538920 927538921 1311738121 } { 5406093003 5406093004 7645370045 } { 31509019100 31509019101 44560482149 } { 183648021599 183648021600 259717522849 } }
</pre>
=={{header|Ruby}}==
<syntaxhighlight lang="ruby">
require 'openssl'
def prime?(n) = OpenSSL::BN.new(n).prime?
pell = Enumerator.new do |y|
pe = [0,1]
loop{y << pe[0]; pe << pe[1]*2 + pe.shift}
end
pell_lucas = Enumerator.new do |y|
pe = [2,2]
loop{y << pe[0]; pe << pe[1]*2 + pe.shift}
end
n = 20
puts "First #{n} Pell numbers: #{pell.first(n).to_a.inspect}"
puts "\nFirst #{n} Pell-Lucas numbers: #{pell_lucas.first(n).to_a.inspect}"
n = 15
sqrt2 = pell.each_cons(2).lazy.map{|p1,p2| Rational(p1+p2, p2)}.take(n).to_a
puts "\nThe first #{n} rational approximations of √2 (#{Math.sqrt(2)}) are:"
sqrt2.each{|n| puts "%-16s ≈ %-18s\n"% [n, n.to_f]}
puts "\nThe first #{n} Pell primes with index are:"
primes = pell.with_index.lazy.select{|n, i|prime?(i) && prime?(n)}.first(n)
primes.each {|prime, i| puts "#{i.to_s.ljust(3)} #{prime}"}
puts "\nThe first #{n} NSW numbers:"
puts pell.first(2*n).each_slice(2).map(&:sum).join(", ")
puts "\nFirst #{n} near isosceles right triangles:"
sum = 1
pell.take(n*2+2).each_slice(2) do |p1,p2|
next if p1 == 0
sum += p1
puts "#{sum}, #{sum+1}, #{p2}"
sum += p2
end
</syntaxhighlight>
{{out}}
<pre>First 20 Pell numbers: [0, 1, 2, 5, 12, 29, 70, 169, 408, 985, 2378, 5741, 13860, 33461, 80782, 195025, 470832, 1136689, 2744210, 6625109]
First 20 Pell-Lucas numbers: [2, 2, 6, 14, 34, 82, 198, 478, 1154, 2786, 6726, 16238, 39202, 94642, 228486, 551614, 1331714, 3215042, 7761798, 18738638]
The first 15 rational approximations of √2 (1.4142135623730951) are:
1/1 ≈ 1.0
3/2 ≈ 1.5
7/5 ≈ 1.4
17/12 ≈ 1.4166666666666667
41/29 ≈ 1.4137931034482758
99/70 ≈ 1.4142857142857144
239/169 ≈ 1.4142011834319526
577/408 ≈ 1.4142156862745099
1393/985 ≈ 1.4142131979695431
3363/2378 ≈ 1.4142136248948696
8119/5741 ≈ 1.4142135516460548
19601/13860 ≈ 1.4142135642135643
47321/33461 ≈ 1.4142135620573204
114243/80782 ≈ 1.4142135624272734
275807/195025 ≈ 1.4142135623637995
The first 15 Pell primes with index are:
2 2
3 5
5 29
11 5741
13 33461
29 44560482149
41 1746860020068409
53 68480406462161287469
59 13558774610046711780701
89 4125636888562548868221559797461449
97 4760981394323203445293052612223893281
101 161733217200188571081311986634082331709
167 2964793555272799671946653940160950323792169332712780937764687561
181 677413820257085084326543915514677342490435733542987756429585398537901
191 4556285254333448771505063529048046595645004014152457191808671945330235841
The first 15 NSW numbers:
1, 7, 41, 239, 1393, 8119, 47321, 275807, 1607521, 9369319, 54608393, 318281039, 1855077841, 10812186007, 63018038201
First 15 near isosceles right triangles:
3, 4, 5
20, 21, 29
119, 120, 169
696, 697, 985
4059, 4060, 5741
23660, 23661, 33461
137903, 137904, 195025
803760, 803761, 1136689
4684659, 4684660, 6625109
27304196, 27304197, 38613965
159140519, 159140520, 225058681
927538920, 927538921, 1311738121
5406093003, 5406093004, 7645370045
31509019100, 31509019101, 44560482149
183648021599, 183648021600, 259717522849
</pre>
=={{header|Scala}}==
<syntaxhighlight lang="scala">val pellNumbers: LazyList[BigInt] =
BigInt("0") #:: BigInt("1") #::
(pellNumbers zip pellNumbers.tail).
map{ case (a,b) => 2*b + a }
val pellLucasNumbers: LazyList[BigInt] =
BigInt("2") #:: BigInt("2") #::
(pellLucasNumbers zip pellLucasNumbers.tail).
map{ case (a,b) => 2*b + a }
val pellPrimes: LazyList[BigInt] =
pellNumbers.tail.tail.
filter{ case p => p.isProbablePrime(16) }
val pellIndexOfPrimes: LazyList[BigInt] =
pellNumbers.tail.tail.
filter{ case p => p.isProbablePrime(16) }.
map{ case p => pellNumbers.indexOf(p) }
val pellNSWnumbers: LazyList[BigInt] =
(pellNumbers.zipWithIndex.collect{ case (p,i) if i%2 == 0 => p}
zip
pellNumbers.zipWithIndex.collect{ case (p,i) if i%2 != 0 => p}).
map{ case (a,b) => a + b }
val pellSqrt2Numerator: LazyList[BigInt] =
BigInt(1) #:: BigInt(3) #::
(pellSqrt2Numerator zip pellSqrt2Numerator.tail).
map{ case (a,b) => 2*b + a }
val pellSqrt2: LazyList[BigDecimal] =
(pellSqrt2Numerator zip pellNumbers.tail).
map{ case (n,d) => BigDecimal(n)/BigDecimal(d) }
val pellSqrt2asString: LazyList[String] =
(pellSqrt2Numerator zip pellNumbers.tail).
map{ case (n,d) => s"$n/$d" }
val pellHypotenuse: LazyList[BigInt] =
pellNumbers.tail.tail.zipWithIndex.collect{ case (p,i) if i%2 != 0 => p }
val pellShortLeg: LazyList[BigInt] =
LazyList.from(3,2).map{ case s => pellNumbers.take(s).sum }
val pellTriple: LazyList[(BigInt,BigInt,BigInt)] =
(pellHypotenuse zip pellShortLeg).
map{ case (h,s) => (s,s+1,h)}
// Output
{ println("7 Tasks")
println("-------")
println(pellNumbers.take(10).mkString("1. Pell Numbers: ", ",", "\n"))
println(pellLucasNumbers.take(10).mkString("2. Pell-Lucas Numbers: ", ",", "\n"))
println((pellSqrt2asString zip pellSqrt2).take(10).
map { case (f, d) => s"$f = $d" }.
mkString("3. Square-root of 2 Approximations: \n\n",
"\n", "\n"))
println(pellPrimes.take(10).mkString("4. Pell Primes: \n\n", "\n", "\n"))
println(pellIndexOfPrimes.take(10).mkString("5. Pell Index of Primes: ", ",", "\n"))
println(pellNSWnumbers.take(10).mkString("6. Newman-Shank-Williams Numbers: \n\n", "\n", "\n"))
println(pellTriple.take(10).mkString("7. Near Right-triangle Triples: \n\n", "\n", "\n"))
}
</syntaxhighlight>
{{out}}
<pre style="height:48ex;overflow:scroll;>7 Tasks
-------
1. Pell Numbers: 0,1,2,5,12,29,70,169,408,985
2. Pell-Lucas Numbers: 2,2,6,14,34,82,198,478,1154,2786
3. Square-root of 2 Approximations:
1/1 = 1
3/2 = 1.5
7/5 = 1.4
17/12 = 1.416666666666666666666666666666667
41/29 = 1.413793103448275862068965517241379
99/70 = 1.414285714285714285714285714285714
239/169 = 1.414201183431952662721893491124260
577/408 = 1.414215686274509803921568627450980
1393/985 = 1.414213197969543147208121827411168
3363/2378 = 1.414213624894869638351555929352397
4. Pell Primes:
2
5
29
5741
33461
44560482149
1746860020068409
68480406462161287469
13558774610046711780701
4125636888562548868221559797461449
5. Pell Index of Primes: 2,3,5,11,13,29,41,53,59,89
6. Newman-Shank-Williams Numbers:
1
7
41
239
1393
8119
47321
275807
1607521
9369319
7. Near Right-triangle Triples:
(3,4,5)
(20,21,29)
(119,120,169)
(696,697,985)
(4059,4060,5741)
(23660,23661,33461)
(137903,137904,195025)
(803760,803761,1136689)
(4684659,4684660,6625109)
(27304196,27304197,38613965)</pre>
=={{header|Wren}}==
Line 802 ⟶ 2,454:
{{libheader|Wren-math}}
{{libheader|Wren-fmt}}
<
import "./math" for Int
import "./fmt" for Fmt
Line 867 ⟶ 2,519:
p1 = p2
i = i + 1
}</
{{out}}
| |||