The Museum of HP Calculators
The Museum of HP Calculators
This program is Copyright © 2006 by Jean-Marc Baillard and is used here by permission.
This program is supplied without representation or warranty of any kind. Jean-Marc Baillard and The Museum of HP Calculators therefore assume no responsibility and shall have no liability, consequential or otherwise, of any kind arising from the use of this program material or any part thereof.
Overview
1°) Legendre Polynomials
2°) Laguerre Polynomials
3°) Generalized Laguerre Polynomials
4°) Hermite Polynomials
5°) Chebyshev Polynomials & the "Connaissanse
des Temps"
a) Chebyshev Polynomials of the first
and second kind
b) The "Connaissance des Temps"
6°) Ultraspherical Polynomials
7°) Jacobi's Polynomials
1°) Legendre Polynomials
Formulae: n.Pn(x)
= (2n-1).x.Pn-1(x) - (n-1).Pn-2(x) ;
P0(x) = 1 ; P1(x) = x
Data Registers: /
Flags: /
Subroutines: /
-Synthetic register M may be replaced by any unused data register.
01 LBL "LEG"
02 X<>Y
03 E3
04 /
05 1
06 +
07 STO M
08 LASTX
09 LBL 01
10 STO N
11 RCL Z
12 *
13 STO O
14 -
15 RCL M
16 INT
17 1/X
18 1
19 -
20 *
21 RCL O
22 +
23 RCL N
24 X<>Y
25 ISG M
26 GTO 01
27 CLA
28 END
( 46 bytes / SIZE 000 )
| STACK | INPUTS | OUTPUTS |
| T | / | x |
| Z | / | x |
| Y | n | Pn-1(x) |
| X | x | Pn (x) |
( 0 < n < 999 )
Example: Calculate P7(4.9)
7 ENTER^
4.9 XEQ "LEG" gives P7(4.9)
= 1698444.019 ( in 4 seconds )
and P6(4.9) = 188641.3852 is
in Y-register.
2°) Laguerre Polynomials
Formulae: Ln(x)
= (2n-1-x).Ln-1(x) - (n-1)2.Ln-2(x)
; L0(x) = 1 ; L1(x) = 1 - x
Data Registers: /
Flags: /
Subroutines: /
-Synthetic register M may be replaced by any unused data register.
01 LBL "LAG"
02 X<>Y
03 .1
04 %
05 1
06 +
07 STO M
08 SIGN
09 LBL 01
10 X<>Y
11 RCL M
12 INT
13 DSE X
14 ""
TEXT0 or another NOP instruction like LBL 02 STO
X ABS ... etc ...
15 X^2
16 *
17 RCL M
18 INT
19 ST+ X
20 DSE X
21 R^
22 ST- Y
23 X<> T
24 ST* Y
25 X<> Z
26 -
27 ISG M
28 GTO 01
29 CLA
30 END
( 51 bytes / SIZE 000 )
| STACK | INPUTS | OUTPUTS |
| T | / | x |
| Z | / | x |
| Y | n | Ln-1(x) |
| X | x | Ln (x) |
( 0 < n < 1000 )
Example: Calculate L7 (3.14)
7 ENTER^
3.14 XEQ "LAG" >>>> L7 (3.14)
= -4932.43995 (
in 5 seconds )
X<>Y >>>> L6
(3.14) = -189.9784735
Note: -Some authors divide Ln (x)
by n! In this case simply add
RCL M 1 - INT FACT /
after line 28
( but Y Z and
T-registers now contain the former Ln-1 (x) )
-With this definition, L7 (3.14)
= -0.978658720
3°) Generalized Laguerre Polynomials
Formulae: L0(a)
(x) = 1 ; L1(a) (x) = a+1-x
; n.Ln(a) (x) = (2.n+a-1-x).Ln-1(a)
(x) - (n+a-1).Ln-2(a) (x)
Data Registers: /
Flags: /
Subroutines: /
01 LBL "LANX"
02 STO M
03 CLX
04 E3
05 /
06 1
07 ST- Z
08 +
09 STO N
10 X<>Y
11 STO O
12 CLST
13 SIGN
14 LBL 01
15 X<>Y
16 CHS
17 RCL N
18 INT
19 RCL O
20 +
21 *
22 RCL N
23 INT
24 ST+ X
25 RCL O
26 +
27 RCL M
28 -
29 R^
30 *
31 +
32 RCL N
33 INT
34 /
35 ISG N
36 GTO 01
37 CLA
38 END
( 61 bytes / SIZE 000 )
| STACK | INPUTS | OUTPUTS |
| Z | a | Ln-1(a)(x) |
| Y | n | Ln-1(a)(x) |
| X | x | Ln(a)(x) |
( 0 < n < 1000 , n integer )
Example:
1.4 ENTER^
7 ENTER^
PI XEQ "LANX" >>>> L7(1.4)(Pi)
= 1.688893513 ( 5 seconds )
X<>Y L6(1.4)(Pi) = 2.271353727
Note: Ln(0)(x) = (1/n!).Ln(x)
where Ln(x) = Laguerre's Polynomial
4°) Hermite Polynomials
Formulae: Hn(x)
= 2x.Hn-1(x) - 2(n-1).Hn-2(x) ; H0(x)
= 1 ; H1(x) = 2x
Data Registers: /
Flags: /
Subroutines: /
-Synthetic register M may be replaced by any unused data register.
01 LBL "HMT"
02 X<>Y
03 .1
04 %
05 1
06 +
07 STO M
08 SIGN
09 LBL 01
10 X<>Y
11 RCL M
12 INT
13 DSE X
14 ""
TEXT0 or another NOP instruction like LBL 02 STO
X ABS ... etc ...
15 *
16 CHS
17 X<>Y
18 ST* Z
19 X<> Z
20 +
21 ST+ X
22 ISG M
23 GTO 01
24 CLA
25 END
( 42 bytes / SIZE 000 )
| STACK | INPUTS | OUTPUTS |
| T | / | x |
| Z | / | x |
| Y | n | Hn-1(x) |
| X | x | Hn (x) |
( 0 < n < 1000 )
Example: Calculate H7 (3.14)
7 ENTER^
3.14 XEQ "HMT" >>>> H7 (3.14)
= 73726.24332 (
in 4 seconds )
X<>Y >>>> H6
(3.14) = 21659.28040
5°) Chebyshev Polynomials & the "Connaissance
des Temps"
a) Chebyshev
Polynomials of the first and second kind
Formulae: Tn(x) = 2x.Tn-1(x) - Tn-2(x) ; T0(x) = 1 ; T1(x) = x ( first kind ) & Un(x) = 2x.Un-1(x) - Un-2(x) ; U0(x) = 1 ; U1(x) = 2x ( second kind )
Data Registers: /
Flag: F02
if F02 is clear, "CHB" calculates the Chebyshev polynomials of the first
kind: Tn(x)
if F02 is set, "CHB" calculates the Chebyshev polynomials
of the second kind: Un(x)
Subroutines: /
-Synthetic register M may be replaced by any unused data register.
01 LBL "CHB"
02 STO M
03 ST+ M
04 FS? 02
05 CLX
06 X<>Y
07 E3
08 /
09 1
10 +
11 STO Z
12 SIGN
13 LBL 01
14 RCL M
15 X<>Y
16 ST* Y
17 X<> Z
18 -
19 ISG Z
20 GTO 01
21 CLA
22 END
( 40 bytes / SIZE 000 )
| STACK | INPUTS | OUTPUTS |
| Y | n | Tn-1(x)orUn-1(x) |
| X | x | Tn (x)orUn(x) |
( 0 < n < 1000 )
Example: Compute T7 (0.314) and U7(0.314)
CF 02
7 ENTER^
0.314 XEQ "CHB" >>>> T7
(0.314) = -0.786900700 and
T6 (0.314) = 0.338782777 in Y-register
( in 2 seconds )
SF 02
7 ENTER^
0.314 R/S
>>>> U7 (0.314) = -0.582815681
and U6 (0.314) = 0.649952293
in Y-register
b) The "Connaissance des Temps"
-The "Connaissance des Temps" is a French ephemeride which contains
coefficients of Chebyshev polynomials
for very accurate coordinates of the Sun , the Moon and the
planets.
-If y(t) is a coordinate of a celestial body at a given instant t (
t belonging to the interval [ t0 ; t0 + DT
] ) , y is given by the formula:
y = a0 + a1.T1(x)
+ ....... + an.Tn(x)
where x = -1 + 2(t-t0)/DT
( -1 <= x <= +1 )
and a0 ; a1 ; ....... ; an
are published in the Connaissance des Temps
Data Registers: only the (n+1) registers
used to store the coefficients a0 ; a1
; ....... ; an
Flags: /
Subroutines: /
-Synthetic registers M N O may be replaced by any unused data registers.
01 LBL "CdT"
02 X<>Y
03 STO M
04 SIGN
05 ST+ M
06 STO N
07 RCL Y
08 STO O
09 RCL IND L
10 LBL 01
11 R^
12 ST+ X
13 RCL N
14 ST* Y
15 X<> O
16 -
17 STO N
18 RCL IND M
19 *
20 +
21 ISG M
22 GTO 01
23 CLA
24 END
( 46 bytes )
| STACK | INPUTS | OUTPUTS |
| T | / | x |
| Z | / | x |
| Y | bbb.eee | x |
| X | x | y(x) |
where bbb.eee is the control number of the registers which contain the (n+1) coefficients a0 ; a1 ; ....... ; an
Example: Compute the radius vector of Saturn for 2000 March 12 at 0h TT ( TT = TAI+32.184s )
-We find the following coefficients in the "Connaissance des Temps 2000" ( page B37 ) valid for 2000/01/00 0h to 2000/12/33 0h ( DT = 368 days )
a0 = 9.14765315 ; a1 = -0.03544281 ; a2 = 0.00109597 ; a3 = 0.00002140 ; a4 = 0.00000039 ; a5 = -0.00000083 ( in AU )
-For instance, we store these 6 numbers in registers R10 to R15
( bbb.eee = 10.015 )
-In this example t-t0 = 72 days therefore x = -1 +2*72/368
= -0.6086956522
10.015
ENTER^
-0.6086956522 XEQ "CdT" >>>>
radius vector of Saturn = 9.16896253 AU ( in 2
seconds )
6°) Ultraspherical Polynomials
Formulae: C0(a) (x) = 1 ; C1(a) (x) = 2.a.x ; (n+1).Cn+1(a) (x) = 2.(n+a).x.Cn(a) (x) - (n+2a-1).Cn-1(a) (x) if a # 0
Cn(0) (x) = (2/n).Tn(x)
Data Registers: /
Flags:
F02 if a = 0
| none if a # 0
Subroutines: CHB - Chebyshev Polynomials
if a =0 | -------------
01 LBL "USP"
02 1
03 R^
04 X#0?
05 GTO 00
06 CF 02
07 R^
08 R^
09 XEQ "CHB"
10 ST+ X
11 R^
12 INT
13 /
14 RTN
15 LBL 00
16 -
17 ST+ X
18 STO O
19 RDN
20 STO M
21 CLX
22 E3
23 /
24 1
25 +
26 STO N
27 CLST
28 SIGN
29 LBL 01
30 X<>Y
31 RCL O
32 RCL N
33 INT
34 -
35 ST* Y
36 X<> L
37 ST+ X
38 RCL O
39 -
40 RCL M
41 *
42 R^
43 *
44 +
45 RCL N
46 INT
47 /
48 ISG N
49 GTO 01
50 CLA
51 END
( 81 bytes / SIZE 000 )
| STACK | INPUTS | OUTPUTS |
| Z | a | / |
| Y | n | Cn-1(a)(x) |
| X | x | Cn(a)(x) |
( 0 < n < 1000 , n integer )
Cn-1(a)(x) is in Y-register only if a # 0
Example: 1.5 ENTER^
7 PI 1/X XEQ "USP" >>>>
C7(1.5)(1/Pi) = -0.989046386
X<>Y C6(1.5)(1/Pi) =
1.768780932
7°) Jacobi's Polynomials
Formulae: P0(a;b) (x) = 1 ; P1(a;b) (x) = (a-b)/2 + x.(a+b+2)/2
2n(n+a+b)(2n+a+b-2).Pn(a;b)
(x) = [ (2n+a+b-1).(a2-b2) + x.(2n+a+b-2)(2n+a+b-1)(2n+a+b)
] Pn-1(a;b) (x)
- 2(n+a-1)(n+b-1)(2n+a+b) Pn-2(a;b) (x)
Data Registers:
R00 thru R05: temp ( R01 = x )
Flags: /
Subroutines: /
01 LBL"JCP"
02 STO 01
03 CLX
04 E3
05 /
06 1
07 STO 04
08 +
09 STO 00
10 RDN
11 STO 03
12 X<>Y
13 STO 02
14 +
15 2
16 +
17 RCL 01
18 *
19 RCL 02
20 RCL 03
21 -
22 +
23 2
24 LBL 01
25 /
26 ISG 00
27 FS? 30
28 GTO 02
29 X<> 04
30 CHS
31 RCL 00
32 INT
33 RCL 02
34 +
35 STO 05
36 1
37 -
38 *
39 RCL 00
40 INT
41 RCL 03
42 +
43 ST+ 05
44 1
45 -
46 *
47 RCL 05
48 ENTER^
49 ST* Z
50 1
51 -
52 STO 05
53 ST* Y
54 LASTX
55 -
56 *
57 RCL 01
58 *
59 RCL 02
60 X^2
61 RCL 03
62 X^2
63 -
64 RCL 05
65 *
66 +
67 RCL 04
68 *
69 2
70 /
71 +
72 RCL 05
73 1
74 -
75 RCL 02
76 RCL 03
77 +
78 RCL 00
79 INT
80 ST+ Y
81 *
82 *
83 GTO 01
84 LBL 02
85 RCL 04
86 X<>Y
87 END
( 105 bytes / SIZE 006 )
| STACK | INPUTS | OUTPUTS |
| T | a | / |
| Z | b | / |
| Y | n | Pn-1(a;b)(x) |
| X | x | Pn(a;b)(x) |
( 0 < n < 1000 , n integer )
Example:
1.4 ENTER^
1.7 ENTER^
7 ENTER^
PI 1/X XEQ "JCP" >>>> P7(1.4;1.7)(1/Pi)
= -0.322234421 ( in 13 seconds )
X<>Y P6(1.4;1.7)(1/Pi) =
0.538220923
Reference:
Abramowitz and Stegun , "Handbook of Mathematical Functions" - Dover
Publications - ISBN 0-486-61272-4
Go back to the HP-41 software library
Go back to the general software library
Go
back to the main exhibit hall