Or then we use
Horner's method to avoid exponentiation, e.g. for \(n = 5\):
\(
\begin{align}
u &= \frac{x - 1}{2} \\
\\
P_5(x) &= 1 + u \cdot \left(30 + u \cdot \left(210 + u \cdot \left(560 + u \cdot \left(630 + u \cdot 252 \right) \right) \right) \right) \\
\end{align}
\)
It's easy to see that:
\(
\begin{align}
\binom{n}{k-1} &= \binom{n}{k} \frac{k}{n-k+1} \\
\\
\binom {n+k-1}{k-1} &= \binom {n+k}{k} \frac{k}{n + k} \\
\end{align}
\)
This leads to a recursive formula to calculate the coefficient from \(k \to k-1\):
\(
\begin{align}
\binom{n}{k-1} \binom {n+k-1}{k-1} = \binom {n}{k}\binom {n+k}{k} \frac{k}{n-k+1} \frac{k}{n + k}
\end{align}
\)
It is used to calculate the next coefficient from the previous one.
Note: All of these coefficients are integers.
Therefore, the order of operation is important.
We also try to keep the intermediate results small.
Python
Code:
from math import comb
def legendre(n, x):
u = (x - 1) / 2
a = comb(2*n, n)
s = 0
for k in range(n, 0, -1):
s = a + u * s
a = a * k // (n + k) * k // (n - k + 1)
return a + u * s
Example
\( P_5(0.3) \)
5
ENTER
0.3
XEQ "LP"
0.34538625
Programs
In case of the
HP-42S and
HP-15C (and other models as well)
RCL-arithmetic could be used.
Also I tried to avoid stack acrobatics to keep the code similar.
Feel free to optimize the program for your favourite model.
HP-42S
Code:
00 { 52-Byte Prgm }
01▸LBL "LP"
02 1
03 -
04 2
05 ÷
06 STO 00
07 R↓
08 STO 01
09 STO 02
10 RCL 01
11 +
12 RCL 01
13 COMB
14 0
15▸LBL 00
16 RCL 00
17 ×
18 X<>Y
19 +
20 LASTX
21 RCL 02
22 ×
23 RCL 01
24 RCL 02
25 +
26 ÷
27 RCL 02
28 ×
29 RCL 01
30 RCL 02
31 -
32 1
33 +
34 ÷
35 X<>Y
36 DSE 02
37 GTO 00
38 RCL 00
39 ×
40 +
41 END
HP-15C
Code:
000 { }
001 { 42 21 11 } f LBL A
002 { 1 } 1
003 { 30 } −
004 { 2 } 2
005 { 10 } ÷
006 { 44 0 } STO 0
007 { 33 } R⬇
008 { 44 1 } STO 1
009 { 44 2 } STO 2
010 { 45 1 } RCL 1
011 { 40 } +
012 { 43 36 } g LSTΧ
013 { 43 40 } g Cy,x
014 { 0 } 0
015 { 42 21 0 } f LBL 0
016 { 45 0 } RCL 0
017 { 20 } ×
018 { 34 } x↔y
019 { 40 } +
020 { 43 36 } g LSTΧ
021 { 45 2 } RCL 2
022 { 20 } ×
023 { 45 1 } RCL 1
024 { 45 2 } RCL 2
025 { 40 } +
026 { 10 } ÷
027 { 45 2 } RCL 2
028 { 20 } ×
029 { 45 1 } RCL 1
030 { 45 2 } RCL 2
031 { 30 } −
032 { 1 } 1
033 { 40 } +
034 { 10 } ÷
035 { 34 } x↔y
036 { 42 5 2 } f DSE 2
037 { 22 0 } GTO 0
038 { 45 0 } RCL 0
039 { 20 } ×
040 { 40 } +
041 { 43 32 } g RTN
Some models use a specific register
I for loops.
HP-11C
Code:
01-42,21,11 # f LBL A
02- 1 # 1
03- 30 # −
04- 2 # 2
05- 10 # ÷
06- 44 0 # STO 0
07- 33 # R⬇
08- 44 1 # STO 1
09- 44 25 # STO I
10- 45 1 # RCL 1
11- 40 # +
12- 43 36 # g LSTΧ
13- 43 1 # g Cy,x
14- 0 # 0
15-42,21, 0 # f LBL 0
16- 45 0 # RCL 0
17- 20 # ×
18- 34 # x↔y
19- 40 # +
20- 43 36 # g LSTΧ
21- 45 25 # RCL I
22- 20 # ×
23- 45 1 # RCL 1
24- 45 25 # RCL I
25- 40 # +
26- 10 # ÷
27- 45 25 # RCL I
28- 20 # ×
29- 45 1 # RCL 1
30- 45 25 # RCL I
31- 30 # −
32- 1 # 1
33- 40 # +
34- 10 # ÷
35- 34 # x↔y
36- 42 5 # f DSE
37- 22 0 # GTO 0
38- 45 0 # RCL 0
39- 20 # ×
40- 40 # +
41- 43 32 # g RTN
For models that don't have a function to compute binomial coefficients, we can use factorial instead:
\(
\begin{align}
\binom{2n}{n} = \frac{(2n)!}{(n!)^2}
\end{align}
\)
HP-67
Code:
001: 31 25 11 # f LBL A
002: 01 # 1
003: 51 # −
004: 02 # 2
005: 81 # ÷
006: 33 00 # STO 0
007: 35 53 # h R⬇
008: 33 01 # STO 1
009: 35 33 # h ST I
010: 34 01 # RCL 1
011: 61 # +
012: 35 81 # h N!
013: 34 01 # RCL 1
014: 35 81 # h N!
015: 32 54 # g x^2
016: 81 # ÷
017: 00 # 0
018: 31 25 00 # f LBL 0
019: 34 00 # RCL 0
020: 71 # ×
021: 35 52 # h x↔y
022: 61 # +
023: 35 82 # h LST Χ
024: 35 34 # h RC I
025: 71 # ×
026: 34 01 # RCL 1
027: 35 34 # h RC I
028: 61 # +
029: 81 # ÷
030: 35 34 # h RC I
031: 71 # ×
032: 34 01 # RCL 1
033: 35 34 # h RC I
034: 51 # −
035: 01 # 1
036: 61 # +
037: 81 # ÷
038: 35 52 # h x↔y
039: 31 33 # f DSZ
040: 22 00 # GTO 0
041: 34 00 # RCL 0
042: 71 # ×
043: 61 # +
044: 35 22 # h RTN