Posts: 653
Threads: 26
Joined: Aug 2010
I just realized that the 34s rounds displayed values in a way that is different from classic HPs (at least the ones I know of). If in FIX mode the value in X is sufficiently low, the display switches to SCI (or ENG) notation. That's fine. But the thresholds are different:
FIX 2 Display
34s HP
80 [1/x] => 0,01 0,01
2 [/] => 6,25 3 0,01
2 [/] => 3,13 3 3,13 3
100,1 [1/x] => 9,99 3 0,01
100,01 [1/x] => 1,00 2 0,01
The 34s switches as soon as the value itself falls below 0,01  while HPs own devices do so as soon as the rounded value passes this limit (i.e. below 0,005). Is this behaviour intentional? At least a subsequent ROUND seems to return the same results in both cases.
While we're at it: The way SLV works still does not match the manual. The latter claims the user interface is the same as on the 15C, but it isn't. On the 34s, the display setting affects the results, while on the 15C (and 34C and 35s and...) the function always returns the best possible result. Walter, you recently asked for an example and I posted one, but there was no reaction.
Dieter
Posts: 4,587
Threads: 105
Joined: Jul 2005
Dieter,
regarding SLV it depends on Pauli: if we want SLV to operate as in the HP15C (as advertized) then apparently the SW has to be adapted. Else we have to specify how it operates ;)
d:)
Posts: 3,229
Threads: 42
Joined: Jul 2006
The FIX to SCI change over is a troublesome part of probably the hairiest piece of code in the device, I'm not surprised something like this did appear. I've long held a suspicion that HP's display code contained similar quirks that became accepted as normal behaviour.
The code in question is lines 945  976 of trunk/display.c The rounding step would have to be performed before the mode switching checks, but it certainly isn't a simple matter of rearranging the code. The question is how important is make a change like this? An awful lot of testing would have to follow  there are many edge and obscure cases the display code must deal with.
As for the SLV command. Making it not depend on display settings is easy. Change the five x[approx]? tests in trunk/xrom/solve.wp34s to x=? ones instead. Is this desirable behaviour or not?
 Pauli
Posts: 1,216
Threads: 75
Joined: Jun 2011
Quote:
As for the SLV command. Making it not depend on display settings is easy. Change the five x[approx]? tests in trunk/xrom/solve.wp34s to x=? ones instead. Is this desirable behaviour or not?
No, at least not for me. I find it ok that the user is able to determine the solver accuracy with the display setting  it's not always necessary to get a solution with full precision, especially as the solving process may not be the fastest for some complicated expressions.
Franz
Posts: 653
Threads: 26
Joined: Aug 2010
Add a simple round command as the final step in your user routine, and you're done. Easy as that.
But there's a much more important point to this: what exactly are we talking about here? Is it the accuracy of the returned root or is it the tolerance of the function result? In other words: Would FIX 2 in a formatdependent SLV routine return a result with at least two valid decimals, or would this return a root for which the function is approximately zero (when rounded to two decimal places)? That's quite a difference.
Dieter
Edited: 30 Dec 2012, 5:27 p.m.
Posts: 1,216
Threads: 75
Joined: Jun 2011
Quote:
Add a simple round command as the final step in your user routine, and you're done. Easy as that.
It seems you didn't understand. Rounding _after_ solving won't make the solver faster at all. Why should I calculate 12 or 16 digits when I only need 3?
Franz
Posts: 653
Threads: 26
Joined: Aug 2010
Quote:
It seems you didn't understand. Rounding _after_ solving won't make the solver faster at all. Why should I calculate 12 or 16 digits when I only need 3?
Oh, Franz.... #)
The rounding is done within the user routine, not after SLV has returned a fullaccuracy result.
Determine the root of x^25 where f(x)...
...is exactly 0 ...is zero if rounded to 3 digits
001 LBL A 001 LBL A
002 x^2 002 x^2
003 5 003 5
004  004 
005 RTN 005 RDP 03
006 RTN
Got the idea? This method was already mentioned in the 34C and 15C manual. And believe me, on those machines it really made sense  they were slow. Really slow.
Dieter
Posts: 1,216
Threads: 75
Joined: Jun 2011
Quote:
Oh, Franz.... #)
Ok, ok, I've overlooked your part with the "user routine". ;)
But this is something you've mentioned in your previous post: does SLV stop when dx~0 or f(x)~0?
Using ROUND in the user routine would influence the accuracy of f(x) (and dx only indirectly), but I guess the WP34s SLV makes its decision (when to terminate) wrt. dxvalue.
Franz
Posts: 653
Threads: 26
Joined: Aug 2010
Quote:
...but I guess the WP34s SLV makes its decision (when to terminate) wrt. dxvalue.
Take a look at the SLV code (in solve.wp34s) and no guessing is required. ;)
Actually there are three tests that check if the function result (!) is approximately zero (applies to f(a), f(b) and f(new_estimate). And finally SLV checks whether the two last approximations agree within 5 ULP.
So, yes, the current SLV implementation seems to do its displaysettingdependent tests based on f(x). Yes, it also checks if the last two estimates are nearly identical (to all digits except the last), but this is not affected by the display setting.
In other words: the rounding done by SLV can just as well be done within the user's function routine. Simply add one of the round commands as the final step.
Dieter
Edited: 30 Dec 2012, 7:17 p.m. after one or more responses were posted
Posts: 3,229
Threads: 42
Joined: Jul 2006
The tests for f(x) being almost zero can be done via a round in the user's code. In fact, they would be better done there because they only do something useful in FIX mode. I'll make this change.
The termination code avoids rounded compares due to some problems there which resulted in bad termination conditions. I don't recollect the specifics anymore.
That just leaves two rounded comparisons. These are used at the start if a singular guess is given  a search interval needs to be created from such. I suspect any interval creation method would be viable here. Will have a think about what might work nicely.
 Pauli
Posts: 556
Threads: 9
Joined: Jul 2007
Is there a chance after that SLV stops giving pole rather than root as an answer in all cases? I haven't checked it lately (are there any changes to SLV?) but I had couple of examples of doing exactly that. I know ... 'everyone is free to develop their own solver' etc., just asking. I'm not capable of doing it.
Cheers,
Posts: 4,587
Threads: 105
Joined: Jul 2005
Is there a chance asking more specific? I.e. unvealing your 'couple of examples'? That would support the quest ... TIA
d:)
Posts: 556
Threads: 9
Joined: Jul 2007
Yes, there is, I'll post them here first thing in the New Year which is now less than 6 hours away. I was actually hoping PD would remember so I asked him.
All the best for the new 2013 to all!
Posts: 3,229
Threads: 42
Joined: Jul 2006
It certainly is possible to make an attempt to detect poles, it just takes time to code and test and a bit of inspiration to escape the influence of the pole.
 Pauli
Posts: 653
Threads: 26
Joined: Aug 2010
Honestly  if we really wanted the 34s to behave exactly as the 15C, this would require the same code of the Solve function. Many examples presented in the 34C and 15C manual cannot be duplicated on the 34s, simply because it uses a different algorithm. Finding poles or local extrema can be done with HP's original SOLVE implementation, but this will not work with SLV on the 34s.
That's why I would suggest to remove the 15C reference from the 34s manual. The only thing both implementations have in common is the stack layout on entry and on exit. And even this is not exactly the same in both cases.
Dieter
Posts: 653
Threads: 26
Joined: Aug 2010
By the way  there also is a way to make SLV return the root (and not f(x)) with a given accuracy. The idea here is to compare the current estimate with the previous one:
001 LBL A
002 FILL // can be omitted as SLV fills the stack anyway
003 X<> 00
004 
005 RDP 02 // previous  current estimate less than 0,005?
006 x=0?
007 RTN // then return zero, make SLV believe a root was found
008 DROP
009 x^2 // usual function code follows here
010 5
011 
012 RTN
Provide two sufficiently different guesses and store the larger one in R00.
2 [ENTER] 3 [STO] 00
Then start SLV as usual. Even with no accuracy limitiations, e.g. in ALL mode, SLV now returns a result within the specified tolerance:
[ALL] 03
[SLV] [A]
=> 2,24000661691
The exact solution is sqrt(5) = 2,236..., so the error in X is less than 0,005  just as specified.
Dieter
Edited: 31 Dec 2012, 9:16 a.m.
Posts: 4,587
Threads: 105
Joined: Jul 2005
Quote:
That's why I would suggest to remove the 15C reference from the 34s manual.
Done. Will be published in due course.
d:)
Posts: 653
Threads: 26
Joined: Aug 2010
I think I'll have to add a big CAVEAT here. The method I suggested works well if the solver produces a strictly convergent series of estimates, i.e. the current estimate is as least as accurate as the difference between the last two:
x_new  x_exact < x_new  x_old
This is not (always) true for the 34s solver, and it does not apply to various other algorithms eiher.
Here's an example. Let's consider the following simple function:
ln(x)  1 = 0
It has a root at x = e. With 2 and 5 as initial guesses The 34s solver produces the following sequence of estimates:
estimate abs. error

2,00000000000 7,18E1
5,00000000000 2,28E 0
3,00465761171 2,86E1 // here the estimate is 0,286 off
2,99797706383 2,80E1 // but the solver adjusts it by only 0,0068
2,49898853191 2,19E1 // which may lead to the conclusion that
2,71681463363 1,47E3 // the estimate is already that exact
2,85739584873 1,39E1
2,71828330193 1,47E6
2,71828215475 3,26E7
2,71754839419 7,33E4
2,71828182846 9,00E15 // here the estimate is almost exact
2,71828199161 1,63E7 //
2,71828182846 0,00E+0 // and here the estimate evnb is dead on
2,71828191003 8,16E8 // but due to roundoff the function result
2,71828186925 4,08E8 // is not exactly zero and the solver again
2,71828184885 2,04E8 // moves away from the exact solution.
2,71828183866 1,02E8 // Starting with an estimate that's 8,16E8 off
2,71828183356 5,10E9 // it now switches to bisection
2,71828183101 2,55E9 // so that it takes more than 20 iterations
2,71828182973 1,27E9 // to find the final digits
2,71828182910 6,37E10
2,71828182878 3,19E10
2,71828182862 1,59E10
2,71828182854 7,97E11
2,71828182850 3,98E11
2,71828182848 1,99E11
2,71828182847 9,96E12
2,71828182846 4,98E12
2,71828182846 2,49E12
2,71828182846 1,25E12
2,71828182846 6,23E13
2,71828182846 3,11E13
2,71828182846 1,55E13
2,71828182846 7,70E14
2,71828182846 3,90E14
2,71828182846 1,90E14
2,71828182846 9,00E15
2,71828182846 5,00E15
2,71828182846 3,00E15
So the current solver eventually comes to a solution, but it may take a while (which is fine if on the other hand a solution is reliably found in more complex cases). Line 3...5 show why the proposed method of comparing the current estimate with the previous one does not work here.
Here is the sequence returned with 1 and 5 as initial guesses:
estimate abs. error

1,00000000000 1,72E+0
5,00000000000 2,28E+0
3,48533973824 7,67E1
2,24266986912 4,76E1
2,68993254763 2,83E2
3,08763614293 3,69E1 // here the solver moves away from the solution
2,71846362445 1,82E4 // ...comes closer
2,71839999019 1,18E4
2,70416626891 1,41E2 // ...moves away again
2,71828182711 1,35E9 // ...comes very close
2,71834090865 5,91E5 // ...but then moves away again
2,71828182846 0,00E+0 // Here x is exact but f(x) is not zero
2,71828182846 2,00E15 // so the solver prefers this solution
Finally I tried the same with the 35s and 34C solvers.
HP 35s HP34C
estimate abs. error estimate abs. error

1,00000000000 1,72E+0 1,000000000 1,72E+0
5,00000000000 2,28E+0 5,000000000 2,28E+0
3,48533973824 7,67E1 3,485339739 7,67E1
2,75051645695 3,22E2 2,442055911 2,76E1
2,71892340384 6,42E4 2,756337820 3,81E2
1,85946170192 8,59E1 2,720245459 1,96E3
2,71814296203 1,39E4 2,718268115 1,37E5
2,71828182743 1,03E9 2,718281833 5,00E9
2,71828182846 0,00E+0 2,718281828 0,00E+0
2,71828182846 0,00E+0 2,718281829 1,00E9
Which leads to the question if there are some deeper insights on how HP Solve (resp. a particular implementation) works in detail. I found an interesting discussion here in the forum back in 2003, cf. this thread.
Dieter
Posts: 4,587
Threads: 105
Joined: Jul 2005
Hallo Dieter,
Quote:
I found an interesting discussion here in the forum back in 2003, cf. this thread.
Thanks for making us aware of that old thread. Interesting read.
d:)
Posts: 653
Threads: 26
Joined: Aug 2010
Yes, the HP solver (here: the version used in the 34C and 15C) is known to produce reliable results, so maybe we can find out more about the way it works and use this for the 34s.
Anyway, I always thought the 34s solver used some kind of implementation of Brent's algorithm. But in this case it would not require 30 to 40 iterations (!) for a simple function like ln(x) = 1. As shown in my previous post, the current method changes to bisection very early, so it converges very slowly.
I tried a simple Brentstyle implementation in VBA for Excel, and here's the sequence of approximations I got:
Estimate abs. error

2,0000 00000 00000 7,18E1
5,0000 00000 00000 2,28E+0
3,0046 57611 71378 2,86E1
2,7574 19085 23714 3,91E2
2,7180 91677 30546 1,90E4
2,7182 83194 09203 1,37E6
2,7182 81828 50681 4,78E11
2,7182 81828 45905 0,00E+0
2,7182 81828 45905
That's 15+ digits after just 9 function calls (including those for the two initial guesses) in a nice convergent series of approximations.
So what's going on in the 34s? Why does the current solver converge so slowly?
Dieter
Posts: 3,229
Threads: 42
Joined: Jul 2006
The 34S solver uses a variety of different methods. It certainly did use Brent at one stage, I'm not how much of that I removed as it evolved. Brent style quadratic interpolation and secant steps are used as is bisection and sometimes Ridder's method (but only after a bisection step). Bisection is used when the root is bracketed and which ever one of the other methods being used appears to produce a divergent x value  this happens quite a lot unfortunately.
I'm not overly attached to the current solver and would certainly consider replacing it with something better.
 Pauli
Posts: 653
Threads: 26
Joined: Aug 2010
As far as I understand there are various methods that make sure that the interval narrows with each iteration step, i.e. the root (if there is one) is found within the interval initially given  provided the two initial guesses bracket the solution.
I just tried a modified regula falsi algorithm (with Anderson/Björck adjustment)  it finds a root within the given interval quite fast and since the calculation is simple and straightforward (as opposed to the Brent method) one or two more iterations do not matter much.
So I see two major points that should be resolved: first, if the initial guesses do not bracket a root, define a method to find such an interval. Second, find a foolproof way to decide whether a) a solution is found and, more important, b) no solution can be found. I mean, something better than "max_iterations exceeded". ;)
Dieter
Posts: 3,229
Threads: 42
Joined: Jul 2006
That max iteration message is a reside of the C implementation  I was extremely memory constrained and couldn't afford space to do much else.
 Pauli
Posts: 3,229
Threads: 42
Joined: Jul 2006
Have any of the 41 machine code gurus figured out how the advantage pack solver works?
If the algorithm used there is understood, it shouldn't be too hard to implement in the 34S.
 Pauli
Posts: 653
Threads: 26
Joined: Aug 2010
While the mcode specialists analyze the 41's advantage ROM code, here's an interesting document I found. It mentions two HP journal articles by W. Kahan and P. McClellan on the HP solvers built into the 34C and 18C. The proposed solver sounds very promising to me. The PDF can be found here.
What do the experts say?
Dieter
Posts: 3,229
Threads: 42
Joined: Jul 2006
Interesting, I hadn't seen that one.
 Pauli
