The Museum of HP Calculators

HP Forum Archive 16

 Root Finding and the 11CMessage #1 Posted by Les Wright on 26 Jan 2007, 4:28 p.m. Has anyone here had occasion to ever work with the adaptation of Newton's Method found at p. 154 of the 11C Owners Manual? I have been intrigued by the routine because it uses a numerical approximation to the derivative so the user doesn't have to provide it. It goes a little overboard in terms of convergence criteria--the routine stops when both the change in the root between iterations and evaluation of the function at the approximated root fall within specified tolerances. I have ported the routine to the HP41 and have found certain circumstances where the numerical approximation to the derivative eventually leads to division by zero error just before convergence since the function value at x and x + delta are essentially the same. This can be rectified by setting the desired tolerance bigger by at least an order of magnitude. That said, the routine is pretty pokey compared the secant method and modified regula falsi method I have mentioned in another thread, since there are two function evaluations for each iteration rather than just one. It is also unnecessarily long and convoluted--there is no need for two convergence criteria or loop counter. I have reduced the size somewhat by cutting back on this. I have inserted a step that shows or prints interim results as desired. If anyone has had occasion to work with this particular routine from our literature, I would like to hear about your experience with it. Les

 Re: Root Finding and the 11CMessage #2 Posted by Namir on 26 Jan 2007, 7:44 p.m.,in response to message #1 by Les Wright Les, I use a practical approximation to the derivative of f(x): f'(x) = (f(x+h) - f(x)) / h where h = 0.01 * (1 + |x|) Newton's algorithm which is: x = x - f(x) / f'(x) becomes: x(new) = x - h * f(x) / (f(x+h) - f(x)) I usually calculate f(x) first and store it since it appeas twice in the above equation. The iterations continue until the absolute value of the refinement in x (that is the absolute value of h * f(x) / (f(x+h) - f(x))) is less than a tolerance value OR the iteration has exceeded the maximum limit. Namir

 Re: Root Finding and the 11CMessage #3 Posted by Les Wright on 27 Jan 2007, 8:24 p.m.,in response to message #2 by Namir Thanks, Namir, and the email that preceded it. I turned out pretty quickly a rough little routine that works fine both for the 41 and, with a tiny bit of modification, the 11C. It is about 40 steps, and with a little more ingenuity I could make it shorter. I find this interesting since the 11C program in the manual that does the exact same thing is 68 steps. And your approximation to the derivative behaves a lot better and division by zero errors are largely avoided. I have put in a flag 10 dependent way to print out interim results, analagous to IG and SV in the PPC ROM, since I like that sort of thing. Thanks for sharing with me. Les

 Re: Root Finding and the 11CMessage #4 Posted by Gerson W. Barbosa on 28 Jan 2007, 9:27 a.m.,in response to message #3 by Les Wright Hello Les, What do you think of a 35-step RPN program written by Newton himself? ;-) If you don't like fiction like I do, proceed directly to page 9 for the HP-10C program listing. It runs on the HP-33C with no modification. There's no need to credit the real author here as you can tell by the link: Regards, Gerson.

 HP-41 Root Finder (was Re: Root Finding and the 11C)Message #5 Posted by Valentin Albillo on 29 Jan 2007, 8:19 a.m.,in response to message #4 by Gerson W. Barbosa Hi, Gerson: Thanks for directing Les (and forum's visitors in general) to my Newton's Anniversary commemorative short tale, I hope anyone bold enough to actually read it will be suitably pleased, as Datafile readers were at the time. Perhaps Les will also find interesting this Root Finder routine for the HP-41C, which I submitted more than 25 years ago to PPC to be considered for inclusion in the then future PPC ROM. Per my stated objectives, it was intended to be: interactive, with full prompting (8 steps), and also callable as a non-prompting subroutine (35 steps, including the final END) very short (8 + 35 steps) and fast, with user-selectable precision by setting the desired display mode. suitably "ruggerized", so that a zero derivative wouldn't cause and error or premature termination but would simply redirect the search elsewhere controlled to avoid infinite looping in case of no roots, returning back to the calling program with a condition set so that the calling program can test for it and proceed as needed (perhaps try another guess or try) fully optimized, including a worthwhile trick or two (such as using two D-R in succession to multiply the current guess by a very small increment (thus actually being a relative increment) in order to compute a fairly good approximation to the derivative. This is just 2 bytes, instead of "{some constant} *", which would be several bytes for the constant, plus one for the multiplication (or else you'd need to previously store the constant in some register, which would be even more expensive) I'm including here a reasonably legible scan with the listing, the description, the instructions, and the examples. Best regards from V. Edited: 29 Jan 2007, 8:51 a.m.

 Re: HP-41 Root Finder (was Re: Root Finding and the 11C)Message #6 Posted by Les Wright on 30 Jan 2007, 4:32 a.m.,in response to message #5 by Valentin Albillo Thank you so much! It took me just a few minutes to transcribe the routine into a textfile and convert it to both RAW ans barcode format. Here is the listing in text form for conveniece: ``` 01 LBL "RF" 02 "NAME?" 03 AON 04 PROMPT 05 AOFF 06 ASTO 00 07 "X0?" 08 PROMPT 09 LBL "RFP" 10 CF 00 11 STO N 12 50 13 STO O 14 LBL 00 15 RCL N 16 1 17 D-R 18 D-R 19 + 20 XEQ IND 00 21 STO M 22 RCL N 23 XEQ IND 00 24 X=0? 25 GTO 01 26 ST- M 27 RCL M 28 X=0? 29 SIGN 30 / 31 D-R 32 D-R 33 ST- N 34 RND 35 X=0? 36 GTO 01 37 DSE O 38 GTO 00 39 SF 00 40 LBL 01 41 RCL N 42 CLA 43 END ``` I really like being able to encourage faster convergence by using a less precise setting of FIX, which will affect the decisions made by the steps "RND X=0?" that crop up a couple of times. This is sometimes necessary with some functions or some initial guess to get the routine to stop, since due to rounding the last (usually 10th) SD of the root or function value at the root wiggles around a bit between iterations. Choosing FIX 7 or FIX 8 or even something less will fix this, and typically seems to result in convergence to the root that gets you as close to a zero function value as you can given the accumulated rounding error. I wrote my own little routine with what Namir provided, but it is still 40 steps or so without the prompting stuff at the beginning. I have made your routine my favourite Newton routine and have enshrined it on its very own magnetic card. In the tradition of the PPC ROM I actually prefer the nonprompting subroutine. You should be proud that RFP is only three steps longer than the SV roujtine that is in the PPC ROM. Like SV for my own use I have added two steps to RFP (not in the listing above) that will show or print interim results if flag 10 is set. Many thanks for directing us to this! Les I have also noticed that RFP is short enough to be adapted to the 33C while leaving enough space for functions to be solved. For that adaptation I would shorten it further by removing the 50 iteration criterion--actually hard to implement on the 33C without the DSE function. With all of the steps saved there you could insert a PSE if desire to watch for convergence oneself. Edited: 30 Jan 2007, 4:33 a.m.

 Re: HP-41 Root Finder (was Re: Root Finding and the 11C)Message #7 Posted by Les Wright on 30 Jan 2007, 6:10 a.m.,in response to message #6 by Les Wright Quote: I have also noticed that RFP is short enough to be adapted to the 33C while leaving enough space for functions to be solved Actually, maybe not. The 33C doesn't have an RND function, which is crucial to the conciseness of the routine.

 Re: HP-41 Root Finder (was Re: Root Finding and the 11C)Message #8 Posted by Valentin Albillo on 30 Jan 2007, 6:19 a.m.,in response to message #6 by Les Wright Hi, Les: Thanks a lot for your kind words and appreciation of my humble (and old!) root-finding routine for the HP-41C. Despite its simplicity, it certainly delivers and I've used it myself as a subroutine called from a number of programs with perfectly satisfactory results. As I see you're keen on numerical analysis topics such as definite integration and root finding, I fully recommend that you don't miss my future articles on advanced techniques for both, belonging to a new series due next in Datafile. Very powerful stuff, as you shall see. P.S.: BTW, this next January-February issue of Datafile does feature no less than three articles by me, introducing two new regular series, among many other worthwhile contents. I'll post a brief abstract as soon as the issue is released. Best regards from V.

 Re: Root Finding and the 11CMessage #9 Posted by Les Wright on 30 Jan 2007, 6:21 a.m.,in response to message #4 by Gerson W. Barbosa Quote: It runs on the HP-33C with no modification. In the absence of the trusty RND function on the 33C, does the conditional test x=0? compare the ACTUAL value of the X-register to zero, or does it compare the value rounded according to the setting of FIX? I ask because with that excellent routine there will be times when the difference computed at step 21 will not be exactly zero, and the routine won't converge according to that criterion. Maybe I should just try it out myself! Les

 Re: Root Finding and the 11CMessage #10 Posted by Valentin Albillo on 30 Jan 2007, 6:30 a.m.,in response to message #9 by Les Wright Hi again, Les: The routine featured in my short story is perfectly good and optimal for the HP-10C only, where it must cope with the lack of RAM, the lack of subroutine capability (as it must call your f(x) twice per iteration), and the lack of sorely needed conditionals, which must be mimicked by using "X=0?" and such. For any other HP calc model not hindered by such appalling restrictions, this routine is sub-optimal and you would do much better with some adaptation of my HP-41C routine instead, for instance. Best regards from V.

 Re: Root Finding and the 11CMessage #11 Posted by Les Wright on 30 Jan 2007, 6:58 a.m.,in response to message #9 by Les Wright I tried it out, and found that, alas, on the 33C, zero means zero internally, not the displayed value. This equation turns up in Acton's Numerical Methods that Usually Work. On the 33C: ENTER, x^2, 1, +, 2, /, 29.672, /, - This fits just barely after the 34 step routine on the 33C, though I guess one could save steps by prestoring the constant 29.672. This is basically a quadratic with two real roots, approx 0.017 and 59.3. If you start with a guess of 31, the routine won't converge but will wiggle back and forth between 59.32714430 and 59.32714429. Start with a better guess of 59.3, and the routine stops at the correct 10-digit root of 59.32714431. If the 33C had a RND function, the routine could be made to stop by setting FIX 7. The rounding error can be fussy. If I start with the ostensibly even better guess of 59.327, I get the same wiggling as I got with a starting guess of 31. But if I choose an even worse guess of 30, it converges and stops at 59.32714431! Weird. (Actually, not weird, as I know it has to do with how the rounding error accumulates according to how the root is approached.) Thanks for sharing these amazing little routines. I note that they are half the length of the routine in the 11C manual, and I think are more efficient. Les