Post Reply 
(71B) Bisection Plus for the HP-71B
06-11-2021, 06:37 PM (This post was last modified: 06-25-2021 04:24 AM by C.Ret.)
Post: #2
RE: (71B) Bisection Plus for the HP-71B
Hi,

I really appreciate this algorithm that stay in the inputted range.
Today, it serve me a bit seeking for the roots of a bunch of polynomials.

I just modify the original program to add two features :
  1. Function definition no more enter into the program listing (easier customization)
  2. A short cut exiting the loop when an exact root (Y=0) is fortunately found (this only happened with test polynomial with integer coefficients).


As useall, since I modified and enhance a few the original code, I post here to share for further user .

5 DIM F$[49]
10 INPUT "F(X)=",F$;F$ @ INPUT "[A,B]=";A,B @ INPUT "Tol=","1E-8";T @ INPUT "MaxIte=","30";M
12 X=B @ Y2=VAL(F$) @ X=A @ Y1=VAL(F$) @ IF Y1*Y2>0 THEN BEEP @ DISP "Fa & Fb = sign !" @ END
14 FOR N=1 TO M @ P=X @ X=(A+B)/2 @ Y=VAL(F$) @ IF ABS(B-A)<T THEN 24
16 IF Y*Y1<0 THEN S=(Y1-Y)/(A-X) @ I=S*A-Y1 ELSE S=(Y2-Y)/(B-X) @ I=S*B-Y2
18 L=X @ Z=Y @ X=I/S @ Y=VAL(F$) @ IF Y=0 OR ABS(P-X)<T THEN 24
20 IF Z*Y<0 THEN A=L @ Y1=Z @ B=X @ Y2=Y ELSE IF Y*Y2<0 THEN A=X @ Y1=Y ELSE B=X @ Y2=Y
22 DISP X @ NEXT N @ BEEP
24 DISP USING '3ZX,"F("K,")="SDE';N,T*(X DIV T),Y @ END


Usage:
  1. Press [RUN]. The program prompts you to enter the function definition. Be sure to enter the function without BASIC syntax error, using X (uppercase variable) as the independent variable. The previous function definition may be clear by pressing the [ON] key. Key in the new function definition or edit the proposed one and press [END LINE] to terminate the input.
  2. The program prompt you for the range [A,B]. Key in the values of A and the value of B separated by a comma. Press [END LINE] to validate.
  3. Key in the value for the tolerance or use the default one.
  4. Key in the maximum iteration count and press [END LINE].

The program displays intermediate refinement for the guess. Eventually, adjust display speed with DELAY command. This have to be tuned before starting the program.
When it converges or after the maximal iteration count, the program displays the root value in the format :
00n F( x.xxxxxxxx ) = ±#E-0##.

Where:
00n indicate iteration count.
x.xxxxxxxx is the root \( x_0 \) rounded to the tolerance.
±#E-0## is the sign and amplitude of the residual \( y=F(x_0)\approx 0 \)

Another root can be seek for the same function or for a new function, just press [ f ][CONT] or [RUN] to process. All the parameters can be kept or modified as wish... no edition of the program is needed.


Step 1: Root of \( f(x)=e^x-3x^2 \) in \( \left [\,-2\,,\,0\,\right ] \)
Code:
                                   [RUN]
F(X)=_                             [f][EXP(] X ) - 3 * X [g] ^ 2 [END LINE]
[A,B]=_                            -2,0[END LINE]
Tol=1E-8                           [→][→][→][4][END LINE]                     ;  1E-8 change to 1E-4
MaxIte=30                          [END LINE]   
-.275321257597    -.432885068285    -.457442411298    -.458918397673
005 F(-.4589)=+2E-006              [f][CONT]
Step 2: Root of \( f(x)=e^x-3x^2 \) in \( \left [\,0\,,\,3\,\right ] \)
Code:

F(X)=EXP(X)-3*X^2                  [END LINE]
[A,B]=_                            0,3[END LINE]
TOL=1E-8                           [END LINE]
MaxIte=30                          [END LINE]
.458952661568   .883433519492   .909669444836   .91000548197 .910004566002
006 F(.91000757)=+2E-011           [f][CONT]
Step 3: Root of \( f(x)=e^x-3x^2 \) in \( \left [\,2\,,\,4 \right ] \)
Code:
F(X)=EXP(X)-3*X^2                  [END LINE]
[A,B]=                             2,4[END LINE]
TOL=1E-8                           [END LINE]
MaxIte=30                          [END LINE]
3.67759480515   3.72781796368   3.73284101745   3.73307361257   3.73307896662   3.73307902828
007 F(3.73307902)=+1E-010          [OFF]

Edited 25.06.2021 for typos and examples reformat
Find all posts by this user
Quote this message in a reply
Post Reply 


Messages In This Thread
RE: (71B) Bisection Plus for the HP-71B - C.Ret - 06-11-2021 06:37 PM



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