# hp41programs

MinMaxParabola Min-Max Polynomial for the HP-41

Overview

1°)  Min-Max Parabola
2°)  Min-Max Polynomial

1°) Min-Max Parabola

-Given a set of n data points ( x1 , y1 ) ..... ( xn , yn )     ( n > 3 )
the min-max parabola:  y = a.x2 + b.x + c  minimizes the maximum "error"     h = supi  | yi - a xi2 -  b xi - c |    over the entire set of points.

-The following program uses the "exchange method":

1°) The min-max parabola is found for the first 4 points.
2°) The greatest error is computed at all data points for this parabola.
3°) One of the first 4 points is exchanged with a data point at which the greatest error occurs ( so that the errors are still of alternating sign for this new quadruple )

-The process is repeated until the greatest error occurs at one of the first 4 points ( step 3 is skipped if it happens the first time )

Data Registers:           •  R00 = n = number of points                         ( Registers R00 thru R2n are to be initialized before executing "MMP2" )

•  R01 = x1   •  R03 = x2    ..........  •  R2n-1 = xn                    ( it's necessary to have x1 < x2 < x3 < x4 )
•  R02 = y1   •  R04 = y2    ..........  •  R2n = yn

Flag:  F10 is cleared at the end
Subroutines: /

-Lines 122-160-175-186-196-206-214  are three-byte GTOs

 01  LBL "MMP2"   02  LBL 01   03  RCL 06    04  RCL 02           05  -   06  RCL 07   07  RCL 03   08  -   09  *   10  RCL 08   11  RCL 04   12  -   13  RCL 05   14  RCL 01   15  -   16  *   17  -   18  RCL 01   19  RCL 03   20  -   21  RCL 07   22  -   23  RCL 05   24  ST+ Y   25  RCL 01   26  -   27  *   28  RCL 07   29  RCL 03   30  -   31  *   32  /   33  STO M    34  RCL 01   35  RCL 05   36  +   37  * 38  RCL 06   39  RCL 02    40  -   41  RCL 05           42  RCL 01   43  -   44  /   45  X<>Y   46  -   47  STO N   48  RCL 03   49  RCL 05   50  +   51  *   52  RCL 03   53  X^2   54  RCL 05   55  X^2   56  +   57  RCL M    58  *   59  +   60  RCL 04   61  RCL 06   62  +   63  X<>Y   64  -   65  2   66  /   67  STO O   68  RCL 01   69  RCL M   70  *   71  RCL N   72  +   73  RCL 01   74  * 75  +   76  RCL 02    77  -   78  CF 10   79  X<0?   80  SF 10   81  ENTER^   82  ABS   83  RCL 00           84  ST+ X   85   E3   86  /   87  8   88  +   89  STO P    90  STO Q   91  LBL 02   92  ISG P   93  X<0?   94  GTO 03   95  CLX   96  RCL IND P   97  RCL M   98  *   99  RCL N 100  + 101  RCL IND P 102  * 103  RCL O 104  + 105  ISG P 106  RCL IND P 107  - 108  ABS 109  X<=Y? 110  GTO 02 111  LASTX 112  X<>Y 113  RCL P  114  STO Q 115  GTO 02 116  LBL 03 117  X<> Q 118  INT 119  STO L 120  8 121  X=Y? 122  GTO 09  123  R^ 124  X<0? 125  FS?C 10 126  FS? 30 127  SF 10 128  LASTX 129  DSE X 130  STO Y 131  RCL IND X 132  RCL 07 133  X 03 150  X<> 05 151  X<> 07 152  X<> IND T 153  STO 01 154  RCL 02         155  X<> 04 156  X<> 06 157  X<> 08 158  X<> IND L 159  STO 02 160  GTO 01  161  LBL 04 162  FS? 10 163  GTO 00 164  X<> 05 165  X<> 03 166  X<> 01 167  X<> IND T 168  STO 07 169  RCL 02 170  X<> IND L 171  X<> 08 172  X<> 06 173  X<> 04 174  STO 02 175  GTO 01 176  LBL 05 177  FC? 10 178  GTO 06 179  LBL 00 180  RCL 07 181  X<> IND T 182  STO 07 183  RCL 08 184  X<> IND L 185  STO 08 186  GTO 01  187  LBL 06 188  FS? 10 189  GTO 00 190  RCL 05 191  X<> IND T 192  STO 05  193  RCL 06 194  X<> IND L 195  STO 06 196  GTO 01  197  LBL 07 198  FS? 10 199  GTO 00 200  LBL 08 201  X<> IND T 202  STO 01 203  RCL 02 204  X<> IND L 205  STO 02         206  GTO 01  207  LBL 00 208  RCL 03 209  X<> IND T 210  STO 03 211  RCL 04 212  X<> IND L 213  STO 04 214  GTO 01  215  LBL 09 216  R^ 217  RCL O 218  RCL N 219  RCL M 220  CF 10 221  CLA 222  END

( 321 bytes / SIZE 2n+1 )

 STACK INPUTS OUTPUTS T / +/-h Z / c Y / b X / a

Example:     With the following set of n = 11 points    ( 11  STO 00 )

 x 0 1 2 3 4 5 6 7 8 9 10 y 1.2 6 9.7 13.7 16.8 20 21.7 24 24.9 25.7 26

-Store these 22 numbers into

R01      R03     ..............................................................................................      R21
R02      R04     ..............................................................................................      R22      respectively

XEQ "MMP2"  yields          a = -0.25          ( in 30 seconds )
RDN          b =  5
RDN          c = 0.975
RDN    +/- h = -0.275

>>>  whence the equation of the min-max parabola is  y = -0.25 x2 + 5 x + 0.975  and the maximum "error" is  0.275

Notes:

-When the program stops, the order of the points will have probably changed.
-The abscissas can be unequally spaced.
-Due to roundoff errors, it's perhaps not impossible that the algorithm oscillates for ever between 2 very close | h |-values
-Add  VIEW Z  after line 120 to check
-If it ever happens, set flag F21 ( the program will stop each time VIEW Z is executed ) and press XEQ 09 to get the results.

2°) Min-Max Polynomial of degree d

-The same "exchange method" may also be used to find the min-max polynomial of degree d
-The data must contain at least (d+2) points

WARNING:

-The set of n data points ( x1 , y1 ) ..... ( xn , yn ) must be stored so that   x1 < x2 < x3 < ............... < xd+2

Data Registers:           •  R00 = n = number of points < 45                        ( Registers R00 thru R2n are to be initialized before executing "MMP" )

•  R01 = x1   •  R03 = x2    ..........  •  R2n-1 = xn     with     x1 < x2 < x3 < ............... < xd+2
•  R02 = y1   •  R04 = y2    ..........  •  R2n = yn

R89 to R99: temp                           R100 to R105+d^2+5.d = the coefficients of the successive linear systems

Flag:  CF01  <-> pivoting
SF 01 <-> no pivoting

Subroutine:  "LS"  ( cf "Linear and non-linear systems for the HP-41" )

-Lines 121-169-214-235-245-273 are three-byte GTOs

 01 LBL "MMP"  02 STO 97  03 STO 99            04 X^2  05 LASTX  06 5  07 *  08 +  09 105  10 +  11 STO 98   12 2  13 ST+ 97  14 RCL 97  15 1  16 %  17 RCL 98  18 +  19 .1  20 %  21 100  22 +  23 STO 95  24 RCL 97  25 ENTER  26 X^2  27 .1  28 %  29 +  30 DSE Y  31 +  32 100.1  33 +  34 STO 94  35 RCL 00  36 STO 96  37 GTO 12  38 LBL 11  39 RCL IND 89  40 X<> IND 92  41 STO IND 89  42 CLX  43 SIGN 44 ST+ 89  45 ST+ 92  46 RCL IND 89  47 X<> IND 92  48 STO IND 89  49 LBL 12  50 RCL 98            51 RCL 97  52 ST+ X  53 2 E-5  54 +  55 LBL 00  56 RCL IND X  57 STO IND Z  58 RDN  59 DSE Y  60 DSE X  61 GTO 00  62 FRC  63 RCL 97  64 ST+ X  65 DSE X  66 +  67 STO N  68 CLX  69 RCL 99  70 STO O  71 LBL 01  72 CLX  73 RCL IND N  74 RCL O  75 Y^X  76 STO IND Y  77 DSE Y  78 DSE N  79 GTO 01  80 CLX  81 RCL N  82 CHS  83 RCL 97  84 DSE X  85 ST+ X  86 + 87 STO N  88 DSE O  89 GTO 01  90 RCL 97            91 STO T  92 SIGN  93 LBL 02  94 STO IND Z  95 DSE Z  96 DSE L  97 GTO 02  98 LBL 03  99 CHS 100 STO IND Z 101 DSE Z 102 DSE T 103 GTO 03 104 RCL 95 105 0 106 XEQ "LS" 107 RCL 96 108 STO 00 109 CLX 110 STO 91 111 STO 92 112 RCL 00 113 ST+ X 114 .1 115 % 116 RCL 97 117 ST+ X 118 + 119 ISG X 120 X=0? 121 GTO 10 122 STO 93 123 STO M 124 SIGN 125 RCL 94 126 INT 127 + 128 RCL 97 129 - 130 RCL IND X 131 STO 90 132 VIEW X 133 LBL 04 134 RCL 94           135 ENTER 136 CLX 137 LBL 05 138 RCL IND M 139 * 140 RCL IND Y 141 + 142 DSE Y 143 GTO 05 144 ISG M 145 RCL IND M 146 - 147 STO O 148 ABS 149 RCL 91 150 ABS 151 X>Y? 152 GTO 04 153 RCL O 154 STO 91 155 RCL M 156 STO 92 157 LBL 04 158 ISG M 159 GTO 04 160 DSE 92 161 CLX 162 RCL 90 163 ABS 164  E-5 165 + 166 RCL 91 167 ABS 168 X<=Y? 169 GTO 10 170 RCL 93 171 INT 172 STO Y 173 2 174 - 175 RCL IND 92 176 ENTER 177 LBL 06 178 CLX 179 RCL IND Z 180 X Z 183 STO T 184 X<> Z 185 LBL 06 186 DSE Z 187 CLX 188 DSE Z 189 GTO 06 190 R^ 191 2 192 - 193 STO N 194 ABS 195 STO 89           196 RCL 94  197 ENTER 198 CLX 199 LBL 07 200 RCL IND 89 201 * 202 RCL IND Y 203 + 204 DSE Y 205 GTO 07 206 ISG 89 207 CLX 208 RCL IND 89 209 - 210 RCL 91 211 * 212 DSE 89 213 X>0? 214 GTO 11 215 1 216 CHS 217 RCL N 218 X#Y? 219 GTO 14 220 RCL 92  221 INT 222 1 223 + 224 STO Y 225 2 226 - 227 LBL 08 228 RCL IND Y 229 X<> IND Y 230 STO IND Z 231 RDN 232 DSE Y 233 DSE X 234 GTO 08 235 GTO 12 236 LBL 14 237 RCL 93           238 INT 239 2 240 - 241 X=Y? 242 GTO 14 243 LASTX 244 ST+ 89 245 GTO 11 246 LBL 14 247 RCL 01 248 X<> IND 92 249 STO 01 250 RCL 92 251 1 252 + 253 RCL 02 254 X<> IND Y 255 STO 02 256 RCL 93 257 INT 258 DSE X 259 .1 260 % 261 3 262 + 263 1 264 LBL 09 265 RCL IND Y 266 X<> IND Y 267 STO IND Z 268 CLX 269 SIGN 270 + 271 ISG Y 272 GTO 09 273 GTO 12 274 LBL 10 275 RCL 94           276 RCL 00 277 ST+ X 278 1 279 + 280 STO O 281 LBL 13 282 RCL IND Y 283 STO IND Y 284 CLX 285 SIGN 286 + 287 DSE Y 288 GTO 13 289 1 290 - 291  E3 292 / 293 RCL O 294 + 295 RCL IND Y 296 ABS 297 X<>Y 298 CLA 299 CLD 300 END

( 486 bytes / SIZE 106+d^2+5.d )

 STACK INPUTS OUTPUTS Y / | h | X d bbb.eee

Where       d       =  degree of the polynomial
bbb.eee  =  control number of the min-max polynomial  ( in fact bbb = 2 n + 1 )
and        | h |     =  supi  | yi -  p( xi ) |

Example:    With the following set of 7 data points, find the min-max polynomial of degree d = 3

 x 1 2 4 7 9 10 12 y 2 5 3 7 6 12 17

n = 7   STO 00   CF 01

1  STO 01   2  STO 03   4  STO 05   7  STO 07  9  STO 09  10  STO 11  12  STO 13
2  STO 02   5  STO 04   3  STO 06   7  STO 08  6  STO 10  12  STO 12  17  STO 14

SIZE 130  ( or greater )

3  XEQ "MMP"   >>>>    The HP-41 displays the successive h-values
( namely 1.32  -1.506666667  -1.705   1.75 )  and eventually

bbb.eee = 15.018                                                                    ---Execution time = 4mn43s---
X<>Y         | h |    = 1.75

-The coefficients of the min-max polynomial of degree 3 are in R15-R16-R17-R18 and we get:

p(x) = 0.033333333 x3 - 0.45 x2 + 2.016666666 x + 1.75 = x3 / 30 - 0.45 x2 + ( 121 / 60 ) x + 1.75

and  | h | = 1.75

Notes:

-Since R89 is used for temporary data storage, the number of points must be < 45

-When the program stops, the order of the points will have probably changed.

-The abscissas can be unequally spaced.

-Due to roundoff errors and to avoid a possible infinite loop, I've added lines 164-165
-So, the iterations will stop if 2 consecutive h-values differ by less than 0.00001
-Change these lines ( or delete them ) if you prefer another value...

-If you set flag F01 ( SF 01 ), the linear systems will be solved without pivoting.
-With the example above, you will get ( almost ) the same results in 4m15s, thus saving 28 seconds.
-However, there is a risk of more roundoff-errors.

-If d = 2,  "MMP2" is much faster and for d = 1, cf "Min-Max Line for the HP-41" ( though "MMP" will work too... )

Reference:

[1]  F. Scheid -  "Numerical Analysis" - Mc Graw Hill   ISBN 07-055197-9