Pell numbers
You are encouraged to solve this task according to the task description, using any language you may know.
Pell numbers are an infinite sequence of integers that comprise the denominators of the closest rational approximations to the square root of 2 but have many other interesting uses and relationships.
The numerators of each term of rational approximations to the square root of 2 may also be derived from Pell numbers, or may be found by taking half of each term of the related sequence: Pell-Lucas or Pell-companion numbers.
The Pell numbers: 0, 1, 2, 5, 12, 29, 70, etc., are defined by the recurrence relation:
P0 = 0;
P1 = 1;
Pn = 2 × Pn-1 + Pn-2;
Or, may also be expressed by the closed form formula:
Pn = ((1 + √2)n - (1 - √2)n) / (2 × √2);
Pell-Lucas or Pell-companion numbers: 2, 2, 6, 14, 34, 82, etc., are defined by a very similar recurrence relation, differing only in the first two terms:
Q0 = 2;
Q1 = 2;
Qn = 2 × Qn-1 + Qn-2;
Or, may also be expressed by the closed form formula:
Qn = (1 + √2)n + (1 - √2)n;
or
Qn = P2n / Pn;
The sequence of rational approximations to the square root of 2 begins:
1/1, 3/2, 7/5, 17/12, 41/29, ...
Starting from n = 1, for each term, the denominator is Pn and the numerator is Qn / 2 or Pn-1 + Pn.
Pell primes are Pell numbers that are prime. Pell prime indices are the indices of the primes in the Pell numbers sequence. Every Pell prime index is prime, though not every prime index corresponds to a prime Pell number.
If you take the sum S of the first 4n + 1 Pell numbers, the sum of the terms P2n and P2n + 1 will form the square root of S.
For instance, the sum of the Pell numbers up to P5; 0 + 1 + 2 + 5 + 12 + 29 == 49, is the square of P2 + P3 == 2 + 5 == 7. The sequence of numbers formed by the sums P2n + P2n + 1 are known as Newman-Shank-Williams numbers or NSW numbers.
Pell numbers may also be used to find Pythagorean triple near isosceles right triangles; right triangles whose legs differ by exactly 1. E.G.: (3,4,5), (20,21,29), (119,120,169), etc.
For n > 0, each right triangle hypotenuse is P2n + 1. The shorter leg length is the sum of the terms up to P2n + 1. The longer leg length is 1 more than that.
- Task
- Find and show at least the first 10 Pell numbers.
- Find and show at least the first 10 Pell-Lucas numbers.
- Use the Pell (and optionally, Pell-Lucas) numbers sequence to find and show at least the first 10 rational approximations to √2 in both rational and decimal representation.
- Find and show at least the first 10 Pell primes.
- Find and show at least the first 10 indices of Pell primes.
- Find and show at least the first 10 Newman-Shank-Williams numbers
- Find and show at least the first 10 Pythagorean triples corresponding to near isosceles right triangles.
- See also
- Wikipedia: Pell number
- OEIS:A000129 - Pell numbers
- OEIS:A002203 - Companion Pell numbers (Pell-Lucas numbers)
- OEIS:A001333 - Numerators of continued fraction convergents to sqrt(2) (Companion Pell numbers / 2)
- OEIS:A086383 - Prime terms in the sequence of Pell numbers
- OEIS:A096650 - Indices of prime Pell numbers
- OEIS:A002315 - NSW numbers (Newman-Shank-Williams numbers)
- Wikipedia: Pythagorean triple
- Pythagorean triples
Common Lisp[edit]
(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) )))
The first 10 Pell numbers are:
(pell-sequence 10)
(0 1 2 5 12 29 70 169 408 985 2378 5741 13860)
The first 10 Pell-Lucas numbers are:
(pell-lucas-sequence 10)
(2 2 6 14 34 82 198 478 1154 2786 6726 16238 39202)
Rational approximation of square root of 2:
(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)
The same in decimal form:
(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)
First 7 Pell primes and, as second value, the corresponding indices:
(pell-primes 7)
(0 1 2 5 29 5741 33461 44560482149) ;
(0 1 3 9 11 27)
Still missing some of the taks
FreeBASIC[edit]
#define isOdd(a) (((a) and 1) <> 0)
Function isPrime(Byval ValorEval As Integer) As Boolean
If ValorEval < 2 Then Return False
If ValorEval Mod 2 = 0 Then Return ValorEval = 2
If ValorEval Mod 3 = 0 Then Return ValorEval = 3
Dim d As Integer = 5
While d * d <= ValorEval
If ValorEval Mod d = 0 Then Return False Else d += 2
If ValorEval Mod d = 0 Then Return False Else d += 4
Wend
Return True
End Function
Dim As Integer n
Dim As Integer p(0 To 40), pl(0 To 40)
p(0)= 0: p(1) = 1
pl(0) = 2: pl(1) = 2
For n = 2 To 40
p(n) = 2 * p(n-1) + p(n-2)
pl(n) = 2 * pl(n-1) + pl(n-2)
Next n
Print "First 20 Pell numbers: "
For n = 0 To 19 : Print p(n); : Next n
Print !"\n\nFirst 20 Pell-Lucas: "
For n = 0 To 19 : Print pl(n); : Next n
Print !"\n\nFirst 20 rational approximations of sqrt(2) (" & Str(Sqr(2)) & "): "
For n = 1 To 20
Dim As Integer j = pl(n)/2, d = p(n)
Print Using " &/& ~= &"; j; d; j/d
Next n
Print !"\nFirst 6 Pell primes: [for the limitations of the FB standard library]"
Dim as Integer pdx = 2
Dim As Byte c = 0
Dim As Ulongint ppdx(1 to 20)
do
If isPrime(p(pdx)) Then
If isPrime(pdx) Then ppdx(c) = pdx : End If
Print p(pdx)
c += 1
End If
pdx += 1
loop until c = 6
Print !"\nIndices of first 6 Pell primes: [for the limitations of the FB standard library]"
For n = 0 To 5 : Print " "; ppdx(n); : Next n
Dim As Ulongint nsw(0 To 20)
For n = 0 To 19
nsw(n) = p(2*n) + p(2*n+1)
Next n
Print !"\n\nFirst 20 Newman-Shank-Williams numbers: "
For n = 0 To 19 : Print " "; nsw(n); : Next n
Print !"\n\nFirst 20 near isosceles right triangles:"
Dim As Integer i0 = 0, i1 = 1, i2, t = 1, i = 2, found = 0
Do While found < 20
i2 = i1*2 + i0
If isOdd(i) Then
Print Using " [&, &, &]"; t; t+1 ; i2
found += 1
End If
t += i2
i0 = i1 : i1 = i2
i += 1
Loop
Sleep
- Output:
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: 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.414213562373095): 1/1 ~= 1 3/2 ~= 1.5 7/5 ~= 1.4 17/12 ~= 1.416666666666667 41/29 ~= 1.413793103448276 99/70 ~= 1.414285714285714 239/169 ~= 1.414201183431953 577/408 ~= 1.41421568627451 1393/985 ~= 1.414213197969543 3363/2378 ~= 1.41421362489487 8119/5741 ~= 1.414213551646055 19601/13860 ~= 1.414213564213564 47321/33461 ~= 1.41421356205732 114243/80782 ~= 1.414213562427273 275807/195025 ~= 1.414213562363799 665857/470832 ~= 1.41421356237469 1607521/1136689 ~= 1.414213562372821 3880899/2744210 ~= 1.414213562373142 9369319/6625109 ~= 1.414213562373087 22619537/15994428 ~= 1.414213562373096 First 6 Pell primes: [for the limitations of the FB standard library] 2 5 29 5741 33461 44560482149 Indices of first 6 Pell primes: [for the limitations of the FB standard library] 2 3 5 11 13 29 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]
Go[edit]
package main
import (
"fmt"
"math/big"
"rcu"
)
func main() {
p := make([]int64, 40)
p[1] = 1
for i := 2; i < 40; i++ {
p[i] = 2*p[i-1] + p[i-2]
}
fmt.Println("The first 20 Pell numbers are:")
fmt.Println(p[0:20])
q := make([]int64, 40)
q[0] = 2
q[1] = 2
for i := 2; i < 40; i++ {
q[i] = 2*q[i-1] + q[i-2]
}
fmt.Println("\nThe first 20 Pell-Lucas numbers are:")
fmt.Println(q[0:20])
fmt.Println("\nThe first 20 rational approximations of √2 (1.4142135623730951) are:")
for i := 1; i <= 20; i++ {
r := big.NewRat(q[i]/2, p[i])
fmt.Printf("%-17s ≈ %-18s\n", r, r.FloatString(16))
}
fmt.Println("\nThe first 15 Pell primes are:")
p0 := big.NewInt(0)
p1 := big.NewInt(1)
p2 := big.NewInt(0)
two := big.NewInt(2)
indices := make([]int, 15)
for index, count := 2, 0; count < 15; index++ {
p2.Mul(p1, two)
p2.Add(p2, p0)
if rcu.IsPrime(index) && p2.ProbablyPrime(15) {
fmt.Println(p2)
indices[count] = index
count++
}
p0.Set(p1)
p1.Set(p2)
}
fmt.Println("\nIndices of the first 15 Pell primes are:")
fmt.Println(indices)
fmt.Println("\nFirst 20 Newman-Shank_Williams numbers:")
nsw := make([]int64, 20)
for n := 0; n < 20; n++ {
nsw[n] = p[2*n] + p[2*n+1]
}
fmt.Println(nsw)
fmt.Println("\nFirst 20 near isosceles right triangles:")
u0 := 0
u1 := 1
sum := 1
for i := 2; i < 43; i++ {
u2 := u1*2 + u0
if i%2 == 1 {
fmt.Printf("(%d, %d, %d)\n", sum, sum+1, u2)
}
sum += u2
u0 = u1
u1 = u2
}
}
- Output:
The first 20 Pell numbers are: [0 1 2 5 12 29 70 169 408 985 2378 5741 13860 33461 80782 195025 470832 1136689 2744210 6625109] The first 20 Pell-Lucas numbers are: [2 2 6 14 34 82 198 478 1154 2786 6726 16238 39202 94642 228486 551614 1331714 3215042 7761798 18738638] The first 20 rational approximations of √2 (1.4142135623730951) are: 1/1 ≈ 1.0000000000000000 3/2 ≈ 1.5000000000000000 7/5 ≈ 1.4000000000000000 17/12 ≈ 1.4166666666666667 41/29 ≈ 1.4137931034482759 99/70 ≈ 1.4142857142857143 239/169 ≈ 1.4142011834319527 577/408 ≈ 1.4142156862745098 1393/985 ≈ 1.4142131979695431 3363/2378 ≈ 1.4142136248948696 8119/5741 ≈ 1.4142135516460547 19601/13860 ≈ 1.4142135642135642 47321/33461 ≈ 1.4142135620573205 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.4142135623730964 The first 15 Pell primes are: 2 5 29 5741 33461 44560482149 1746860020068409 68480406462161287469 13558774610046711780701 4125636888562548868221559797461449 4760981394323203445293052612223893281 161733217200188571081311986634082331709 2964793555272799671946653940160950323792169332712780937764687561 677413820257085084326543915514677342490435733542987756429585398537901 4556285254333448771505063529048046595645004014152457191808671945330235841 Indices of the first 15 Pell primes are: [2 3 5 11 13 29 41 53 59 89 97 101 167 181 191] 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)
Haskell[edit]
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
- Output:
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)
J[edit]
As detailed in the task description, there's a variety of ways to compute these values.
For example:
nextPell=: , 1 2+/ .*_2&{. NB. pell, list extender
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:
nextPell^:9(0 1)
0 1 2 5 12 29 70 169 408 985 2378
Pn i.11
0 1 2 5 12 29 70 169 408 985 2378
nextPell^:9(2 2)
2 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 i.11
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 P2n/Pn) doesn't get the first element of the pell lucas sequence right. We could fix this by changing the definition:
QN=: 2 >. +: %&Pn ]
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}. +//.@,:~) nextPell^:9(0 1)
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:
10{.(#~ 1&p:)nextPell^:99(0 1x)
2 5 29 5741 33461 44560482149 1746860020068409 68480406462161287469 13558774610046711780701 4125636888562548868221559797461449
Their indices are:
10{.I. 1&p:nextPell^:99(0 1x)
2 3 5 11 13 29 41 53 59 89
The NSW numbers are the sums of (non-overlapping) pairs of pell numbers, or:
_2 +/\ nextPell^:20(0 1x)
1 7 41 239 1393 8119 47321 275807 1607521 9369319 54608393
The first ten pell based pythogorean triples would be:
}.(21$1 0)#|:(}.,~0 1+/+/\@}:)nextPell^:(20)0 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
Julia[edit]
using Primes
function pellnumbers(wanted)
pells = [0, 1]
wanted < 3 && return pells[1:wanted]
while length(pells) < wanted
push!(pells, 2 * pells[end] + pells[end - 1])
end
return pells
end
function pelllucasnumbers(wanted)
pelllucas = [2, 2]
wanted < 3 && return pelllucas[1:wanted]
while length(pelllucas) < wanted
push!(pelllucas, 2 * pelllucas[end] + pelllucas[end - 1])
end
return pelllucas
end
function pellprimes(wanted)
i, lastpell, lastlastpell, primeindices, pellprimes = 1, big"1", big"0", Int[], BigInt[]
while length(primeindices) < wanted
pell = 2 * lastpell + lastlastpell
i += 1
if isprime(pell)
push!(primeindices, i)
push!(pellprimes, pell)
end
lastpell, lastlastpell = pell, lastpell
end
return primeindices, pellprimes
end
function approximationsqrt2(wanted)
p, q = pellnumbers(wanted + 1), pelllucasnumbers(wanted + 1)
return map(r -> "$r ≈ $(Float64(r))", [(q[n] // 2) // p[n] for n in 2:wanted+1])
end
function newmanshankwilliams(wanted)
pells = pellnumbers(wanted * 2 + 1)
return [pells[2i - 1] + pells[2i] for i in 1:wanted]
end
function nearisosceles(wanted)
pells = pellnumbers((wanted + 1) * 2 + 1)
return map(x -> (last(x), last(x) + 1, first(x)),
[(pells[2i], sum(pells[1:2i-1])) for i in 2:wanted+1])
end
function printrows(title, vec, columnsize = 8, columns = 10, rjust=false)
println(title)
for (i, n) in enumerate(vec)
print((rjust ? lpad : rpad)(n, columnsize), i % columns == 0 ? "\n" : "")
end
println()
end
printrows("Twenty Pell numbers:", pellnumbers(20))
printrows("Twenty Pell-Lucas numbers:", pelllucasnumbers(20))
printrows("Twenty approximations of √2:", approximationsqrt2(20), 44, 2)
pindices, pprimes = pellprimes(15)
printrows("Fifteen Pell primes:", pprimes, 90, 1)
printrows("Fifteen Pell prime zero-based indices:", pindices, 4, 15)
printrows("Twenty Newman-Shank-Williams numbers:", newmanshankwilliams(20), 17, 5)
printrows("Twenty near isosceles triangle triplets:", nearisosceles(20), 52, 2)
- Output:
Twenty Pell numbers: 0 1 2 5 12 29 70 169 408 985 2378 5741 13860 33461 80782 195025 470832 1136689 2744210 6625109 Twenty Pell-Lucas numbers: 2 2 6 14 34 82 198 478 1154 2786 6726 16238 39202 94642 228486 551614 1331714 3215042 7761798 18738638 Twenty approximations of √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 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.414213562373142 9369319//6625109 ≈ 1.414213562373087 22619537//15994428 ≈ 1.4142135623730965 Fifteen Pell primes: 2 5 29 5741 33461 44560482149 1746860020068409 68480406462161287469 13558774610046711780701 4125636888562548868221559797461449 4760981394323203445293052612223893281 161733217200188571081311986634082331709 2964793555272799671946653940160950323792169332712780937764687561 677413820257085084326543915514677342490435733542987756429585398537901 4556285254333448771505063529048046595645004014152457191808671945330235841 Fifteen Pell prime zero-based indices: 2 3 5 11 13 29 41 53 59 89 97 101 167 181 191 Twenty Newman-Shank-Williams numbers: 1 7 41 239 1393 8119 47321 275807 1607521 9369319 54608393 318281039 1855077841 10812186007 63018038201 367296043199 2140758220993 12477253282759 72722761475561 423859315570607 Twenty near isosceles triangle triplets: (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)
Mathematica/Wolfram Language[edit]
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]
- Output:
{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}}
Nim[edit]
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
- Output:
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)
Perl[edit]
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;
- Output:
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)
Phix[edit]
with javascript_semantics sequence p = {0,1}, pl = {2,2} for i=2 to 41 do p &= 2*p[i]+p[i-1] pl &= 2*pl[i]+pl[i-1] end for printf(1,"First 20 Pell numbers: %s\n",{join_by(p[1..20],1,20," ",fmt:="%d")}) printf(1,"First 20 Pell-Lucas: %s\n",{join_by(pl[1..20],1,20," ",fmt:="%d")}) printf(1,"First 20 rational approximations of sqrt(2) (%.16f):\n",{sqrt(2)}) for i=2 to 21 do integer n = pl[i]/2, d = p[i] printf(1,"%d/%d ~= %.16g\n", {n,d,n/d}) end for printf(1,"\nFirst 20 Pell primes:\n") include mpfr.e mpz {p0,p1,p2} = mpz_inits(3,{0,1,0}) sequence ppdx = {} integer pdx = 2 while length(ppdx)<20 do mpz_mul_si(p2,p1,2) mpz_add(p2,p2,p0) if is_prime(pdx) and mpz_prime(p2) then printf(1,"%s\n",mpz_get_short_str(p2)) ppdx = append(ppdx,sprintf("%d",pdx)) end if pdx += 1 mpz_set(p0,p1) mpz_set(p1,p2) end while printf(1,"\nIndices of first 20 Pell primes: %s\n",join(ppdx," ")) sequence nsw = {} for n=1 to 20 do nsw = append(nsw,sprintf("%d",p[2*n]+p[2*n-1])) end for nsw[8..-3] = {"..."} printf(1,"\nFirst 20 Newman-Shank-Williams numbers: %s\n",{join(nsw," ")}) printf(1,"\nFirst 20 near isosceles right triangles:\n") for i=4 to 42 by 2 do atom side = sum(p[1..i-1]), hypot = p[i] printf(1,"[%d, %d, %d]\n", {side,side+1,hypot}) end for
- Output:
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: 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.4142135623730951): 1/1 ~= 1 3/2 ~= 1.5 7/5 ~= 1.4 17/12 ~= 1.416666666666667 41/29 ~= 1.413793103448276 99/70 ~= 1.414285714285714 239/169 ~= 1.414201183431953 577/408 ~= 1.41421568627451 1393/985 ~= 1.414213197969543 3363/2378 ~= 1.41421362489487 8119/5741 ~= 1.414213551646055 19601/13860 ~= 1.414213564213564 47321/33461 ~= 1.41421356205732 114243/80782 ~= 1.414213562427273 275807/195025 ~= 1.414213562363799 665857/470832 ~= 1.41421356237469 1607521/1136689 ~= 1.414213562372821 3880899/2744210 ~= 1.414213562373142 9369319/6625109 ~= 1.414213562373087 22619537/15994428 ~= 1.414213562373097 First 20 Pell primes: 2 5 29 5741 33461 44560482149 1746860020068409 68480406462161287469 13558774610046711780701 4125636888562548868221559797461449 4760981394323203445293052612223893281 161733217200188571081311986634082331709 29647935552727996719...32712780937764687561 (64 digits) 67741382025708508432...87756429585398537901 (69 digits) 45562852543334487715...91808671945330235841 (73 digits) 54971607658948646301...49971207793263738989 (200 digits) 14030291214037674827...9886357032327271649 (356 digits) 24348043146521993819...3590369653438042689 (466 digits) 34643489561492944482...62121115635939797709 (498 digits) 32074710952523740376...14678081652481281009 (521 digits) 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 ... 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]
Quackery[edit]
isprime
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.
[ $ "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 ]
- Output:
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 ]
Raku[edit]
my $pell = cache lazy 0, 1, * + * × 2 … *;
my $pell-lucas = lazy 2, 2, * + * × 2 … *;
my $upto = 20;
say "First $upto Pell numbers:\n" ~ $pell[^$upto];
say "\nFirst $upto Pell-Lucas numbers:\n" ~ $pell-lucas[^$upto];
say "\nFirst $upto rational approximations of √2 ({sqrt(2)}):\n" ~
(1..$upto).map({ sprintf "%d/%d - %1.16f", $pell[$_-1] + $pell[$_], $pell[$_], ($pell[$_-1]+$pell[$_])/$pell[$_] }).join: "\n";
say "\nFirst $upto Pell primes:\n" ~ $pell.grep(&is-prime)[^$upto].join: "\n";
say "\nIndices of first $upto Pell primes:\n" ~ (^∞).grep({$pell[$_].is-prime})[^$upto];
say "\nFirst $upto Newman-Shank-Williams numbers:\n" ~ (^$upto).map({ $pell[2 × $_, 2 × $_+1].sum });
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;
- Output:
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.4142135623730951): 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)
RPL[edit]
Let's start first with a universal Pell sequence generator, which takes a smaller sequence as input:
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
- Output:
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 }
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
- Output:
{ '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 }
To look for prime Pell numbers, we need a fast Pell number generator and also a prime number checker, available 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
TASK4
returns a list of complex numbers, the real part being the prime number and the imaginary part being the indice.
- Output:
1: { (2,2) (5,3) (29,5) (5741,11) (33461,13) (44560482149,29) }
Unfortunately, the emulator's watchdog timer has refused to let us look beyond P29.
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
- Output:
1: { 1 7 41 239 1393 8119 47321 275807 1607521 9369319 54608393 318281039 1855077841 10812186007 63018038201 }
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 ∑LIST
, which sums the items of a list, and replacing LASTARG
with LAST
.
- Output:
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 } }
Scala[edit]
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"))
}
- Output:
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)
Wren[edit]
import "./big" for BigInt, BigRat
import "./math" for Int
import "./fmt" for Fmt
var p = List.filled(40, 0)
p[0] = 0
p[1] = 1
for (i in 2..39) p[i] = 2 * p[i-1] + p[i-2]
System.print("The first 20 Pell numbers are:")
System.print(p[0..19].join(" "))
var q = List.filled(40, 0)
q[0] = 2
q[1] = 2
for (i in 2..39) q[i] = 2 * q[i-1] + q[i-2]
System.print("\nThe first 20 Pell-Lucas numbers are:")
System.print(q[0..19].join(" "))
System.print("\nThe first 20 rational approximations of √2 (1.4142135623730951) are:")
for (i in 1..20) {
var r = BigRat.new(q[i]/2, p[i])
Fmt.print("$-17s ≈ $-18s", r, r.toDecimal(16, true, true))
}
System.print("\nThe first 15 Pell primes are:")
var p0 = BigInt.zero
var p1 = BigInt.one
var indices = List.filled(15, 0)
var count = 0
var index = 2
var p2
while (count < 15) {
p2 = p1 * BigInt.two + p0
if (Int.isPrime(index) && p2.isProbablePrime(10)) {
System.print(p2)
indices[count] = index
count = count + 1
}
index = index + 1
p0 = p1
p1 = p2
}
System.print("\nIndices of the first 15 Pell primes are:")
System.print(indices.join(" "))
System.print("\nFirst 20 Newman-Shank_Williams numbers:")
var nsw = List.filled(20, 0)
for (n in 0..19) nsw[n] = p[2*n] + p[2*n+1]
Fmt.print("$d", nsw)
System.print("\nFirst 20 near isosceles right triangles:")
p0 = 0
p1 = 1
var sum = 1
var i = 2
while (i < 43) {
p2 = p1 * 2 + p0
if (i % 2 == 1) {
Fmt.print("($d, $d, $d)", sum, sum + 1, p2)
}
sum = sum + p2
p0 = p1
p1 = p2
i = i + 1
}
- Output:
The first 20 Pell numbers are: 0 1 2 5 12 29 70 169 408 985 2378 5741 13860 33461 80782 195025 470832 1136689 2744210 6625109 The first 20 Pell-Lucas numbers are: 2 2 6 14 34 82 198 478 1154 2786 6726 16238 39202 94642 228486 551614 1331714 3215042 7761798 18738638 The first 20 rational approximations of √2 (1.4142135623730951) are: 1/1 ≈ 1.0000000000000000 3/2 ≈ 1.5000000000000000 7/5 ≈ 1.4000000000000000 17/12 ≈ 1.4166666666666667 41/29 ≈ 1.4137931034482759 99/70 ≈ 1.4142857142857143 239/169 ≈ 1.4142011834319527 577/408 ≈ 1.4142156862745098 1393/985 ≈ 1.4142131979695431 3363/2378 ≈ 1.4142136248948696 8119/5741 ≈ 1.4142135516460547 19601/13860 ≈ 1.4142135642135642 47321/33461 ≈ 1.4142135620573205 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.4142135623730964 The first 15 Pell primes are: 2 5 29 5741 33461 44560482149 1746860020068409 68480406462161287469 13558774610046711780701 4125636888562548868221559797461449 4760981394323203445293052612223893281 161733217200188571081311986634082331709 2964793555272799671946653940160950323792169332712780937764687561 677413820257085084326543915514677342490435733542987756429585398537901 4556285254333448771505063529048046595645004014152457191808671945330235841 Indices of the first 15 Pell primes are: 2 3 5 11 13 29 41 53 59 89 97 101 167 181 191 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)