(42S) Newton Polynomial Interpolation
03-09-2019, 04:49 PM
Post: #1
 Thomas Klemm Senior Member Posts: 1,447 Joined: Dec 2013
(42S) Newton Polynomial Interpolation
This program uses divided differences to determine the Newton polynomial:

$$N(x)=[y_{0}]+[y_{0},y_{1}](x-x_{0})+\cdots +[y_{0},\ldots ,y_{k}](x-x_{0})(x-x_{1})\cdots (x-x_{{k-1}})$$

Definition of the Polynomial

The polynomial of degree $$n$$ is defined by $$n+1$$ points $$P_i=(x_i, y_i)$$ for $$i \in \{0, \cdots, n\}$$, where all $$x_i$$ are different.

This program can be used to enter the points $$P_i$$:
Code:
00 { 58-Byte Prgm } 01▸LBL "P+" 02 RCL 00 03 1 04 + 05 STO 00 06 3.002 07 × 08 STO 01 09 R↓ 10▸LBL 00 11 STO IND 01 12 DSE 01 13 GTO 01 14 R↓ 15 STO IND 01 16 RTN 17▸LBL 01 18 RCL 01 19 2 20 - 21 R↓ 22 RCL- IND ST T 23 RCL 00 24 2 25 × 26 RCL- 01 27 2 28 × 29 X<> ST T 30 RCL- IND ST T 31 ÷ 32 GTO 00 33 END

Initialisation

Before you enter points, make sure the counter is reset:

0
STO 00

Example:

Find a quadratic polynomial given these 3 points: $$P_0(-5, 12)$$, $$P_1(1, 13)$$ and $$P_2(2, 11)$$

0 STO 00

-5 ENTER 12
XEQ "P+"

1 ENTER 13
XEQ "P+"

2 ENTER 11
XEQ "P+"

These are the coefficients of the Newton polynomial:

R03: 12.000000
R05:  0.166667
R07: -0.309524

$$f(x) = 12 + (x+5)(\frac{1}{6} - (x-1)\frac{13}{42})$$

Interpolation of the Polynomial

This program can be used to evaluate the polynomial:
Code:
00 { 37-Byte Prgm } 01▸LBL "f(x)" 02 1.001 03 RCL 00 04 2 05 × 06 + 07 STO 01 08 CLX 09▸LBL 00 10 RCL IND 01 11 X<>Y 12 DSE 01 13 R↑ 14 RCL- IND 01 15 × 16 + 17 DSE 01 18 GTO 00 19 END

Example:

Evaluate the polynomial at $$x=0.5$$.

0.5
XEQ "f(x)"

13.7679

An advantage of this method is that the polynomial can be evaluated before all points have been entered.
So you can add more points if the interpolation is not considered sufficient.

Example:

Interpolation of $$\sin(x)$$ using values at 0°, 45° and 90°.

0 STO 00

0 ENTER
XEQ "P+"

45 ENTER 0.5 √x
XEQ "P+"

90 ENTER 1
XEQ "P+"

Compare the interpolation at 23° with the correct value:

23
XEQ "f(x)"
0.4132

23
SIN
0.3907

%CH
-5.4289

30 ENTER 0.5
XEQ "P+"

23
XEQ "f(x)"
0.3913

23
SIN
%CH
-0.1397

60 ENTER 0.75 √x
XEQ "P+"

23
XEQ "f(x)"
0.3906

23
SIN
%CH
0.0224

Registers

This is the situation with 3 points:

R00: 3 … number of points
R01: k … index
R02: $$x_0$$
R03: $$[y_0]$$
R04: $$x_1$$
R05: $$[y_0, y_1]$$
R06: $$x_2$$
R07: $$[y_0, y_1, y_2]$$
R08: $$[y_1, y_2]$$
R09: $$[y_2]$$

When the next point is entered, some of the values at the end are overwritten:

R00: 4 … number of points
R01: k … index
R02: $$x_0$$
R03: $$[y_0]$$
R04: $$x_1$$
R05: $$[y_0, y_1]$$
R06: $$x_2$$
R07: $$[y_0, y_1, y_2]$$
R08: $$x_3$$
R09: $$[y_0, y_1, y_2, y_3]$$
R10: $$[y_1, y_2, y_3]$$
R11: $$[y_2, y_3]$$
R12: $$[y_3]$$

Cheers
Thomas

PS: For 3 points I've already posted programs for the HP-67 and HP-15C.
 « Next Oldest | Next Newest »

User(s) browsing this thread: 1 Guest(s)