Hermite Formula

# Osculating Polynomial for the HP-41

Overview

-An osculating polynomial p(x) takes the same values yi  as a given function f at n specified arguments xi     p(xi) = yi = f(xi)   for   i = 1,2,....,n
but its first derivative also verifies:   p'(xi) = y'i = f '(xi)   for   i = 1,2,....,n
-The following program employs Hermite's formula:

p(x) = SUMi=1,2,....,n [ 1 - 2.L'i(xi).(x-xi) ].Li2(x).yi  +  SUMi=1,2,....,n (x-xi).Li2(x).y'i

where    Li(x) = ( x - x1 ) ...... ( x - xi-1 ) ( x - xi+1 ) ...... ( x - xn ) / [ ( xi - x1 ) ...... ( xi - xi-1 ) ( xi - xi+1 ) ( xi - xn ) ]   are the Lagrange's multiplier functions.

-This is the unique polynomial of degree < 2n  that verifies  p(xi) = yi  and  p'(xi) = y'i   for   i = 1,2,....,n
-So, we can hope a better interpolation than with the simple collocation ( Lagrange ) polynomial.

Program Listing

Data Registers:         R00 = bbb.eee                                                           ( Registers Rbbb thru Reee are to be initialized before executing "HPI"  )

•  Rbbb   = x1   •  Rbb+3 = x2  ..........   •  Ree-2 = xn
•  Rbb+1 = y1   •  Rbb+4 = y2  ..........   •  Ree-1 = yn
•  Rbb+2 = y'1  •  Rbb+5 = y'2  ..........   •  Reee  = y'n
Flags: /
Subroutines: /

 01  LBL "HPI"     02  CLA 03  STO M 04  X<>Y 05  STO 00 06  STO O 07  LBL 01 08  CLX 09  STO Q 10  RCL 00 11  STO N 12  SIGN 13  LBL 02 14  RCL M 15  RCL IND O 16  RCL IND N 17  ST- Z 18  - 19  X#0? 20  GTO 03 21  X<> Z 22  GTO 04 23  LBL 03 24  1/X 25  ST- Q 26  * 27  * 28  LBL 04 29  ISG N 30  ISG N 31  ISG N 32  GTO 02 33  X^2 34  STO Y 35  RCL M 36  RCL IND O 37  - 38  ST* Y 39  RCL Q 40  ST+ X 41  * 42  R^ 43  ST* Y 44  + 45  ISG O 46  RCL IND O 47  * 48  X<>Y 49  ISG O 50  RCL IND O 51  * 52  + 53  ST+ P 54  ISG O 55  GTO 01 56  RCL 00 57  RCL M         58  SIGN 59  X<> P 60  CLA 61  END

( 101 bytes )

 STACK INPUTS OUTPUTS Y bbb.eee bbb.eee X x p(x) L / x

where  bbb.eee is the control number of the array  ( bbb > 000 )

Example:

Evaluate  p(6) & p(8)  for the osculating polynomial of degree < 10  that takes the values prescribed below:

 xi 1 2 4 7 10 yi 1 4 6 7 5 y'i 3 2 1 -1 -2

For instance, you store these 15 numbers  into

R01      R04     R07      R10     R13
R02      R05     R08      R11     R14
R03      R06     R09      R12     R15          respectively, then

1.015  ENTER^
6     XEQ "HPI"  >>>>    p(6) = 7.505337940    ( in 20 seconds )

RDN   8   R/S  >>>>  p(8) = 5.746750036

Notes:

-If you don't want to use synthetic registers, replace M N O P Q by R01  R02  R03  R04  R05  and choose  bbb > 005
-On the other hand, you can replace R00 by synthetic register a  but in this case, the program must not be called as more than a first level subroutine.

Reference:

  Francis Scheid   "Numerical Analysis"       McGraw-Hill   ISBN 07-055197-9