FYI: Update to my SOLVE/INTEG article (#556)



#2

All --

For those who are interested, I've made some corrections to my article (#556) in the MoHPC Articles Forum. The changes pertain to INTEG input-function uncertainty and estimated error for the pre-Saturn models. This particular material about INTEG (for all RPN-based models), as well as a list of primary references, was added in early August.

I'll be able to answer any questions only through Saturday; then it'll have to wait until the week after next.

http://www.hpmuseum.org/cgi-sys/cgiwrap/hpmuseum/articles.cgi?read=556

-- KS

Edited: 15 Sept 2006, 12:03 a.m.


#3

I know this question gets asked periodically, but I have yet to see the definitive answer so I'll risk the wrath of the long-timers and ask it again.

Does anybody know the precise algorithm used by the internal integrator on the 34C and 15C (and/or later models)?

I know it's a Romberg with a non-uniform transformation. That much is documented in Kahan's 1980 article and expanded upon in the PPC ROM manual. However, the PPC ROM routine (which was intended to be similar but not identical to the HP rountine) does not give the same results as the 34C/15C or the 41C AdvantagePac (beyond what can be accounted for by the difference in precision). Hugh posted a routine which apparently gives the same results as the 34C/15C with seven function evaluations, but differs with fifteen or more.

Also, the 34C/15C can produce different values for the same number of function evaluations. I don't mean the 34C and the 15C differ, I mean on the same machine you can get two different values for the same set of function evaluations. Try integrating exp(-x^2) from 0 to 1 in FIX 3 and FIX 4.

Any help would be appreciated.


#4

Hi, Koyoshi --

You posted:

Quote:
Also, the 34C/15C can produce different values for the same number of function evaluations. I don't mean the 34C and the 15C differ, I mean on the same machine you can get two different values for the same set of function evaluations. Try integrating exp(-x^2) from 0 to 1 in FIX 3 and FIX 4.

I got the following results:

FIX     Answer     Estimated error  Elapsed Time     
3 0.746825137 0.0005 23.32 sec
4 0.746823672 0.00005 24.96 sec
9 0.746824133 0.0000000005 88 sec

How do you know that the sets of samples were identical? There was at least a half-second difference in time of execution, although it's difficult to ascribe that to more samples taken.

The estimated errors are as indicated by the formula given in the HP-34C and HP-15C manuals, which I quoted in my article. Only five digits (including the leading zero) are needed to display the value of the integrand under FIX 4, so no "roundoff" was performed.

-- KS


#5

Quote:
How do you know that the sets of samples were identical? There was at least a half-second difference in time of execution, although it's difficult to ascribe that to more samples taken.

I put code like [1] [STO + 0] [Rdown] inside the function being integrated. This gives me the number of function evaluations (assuming you know what was in R0 at the beginning). The set of points at which the function is evaluated is known from the non-uniform transformation. For a given number of evals, the set of points is always the same (before the transformation).

I've noticed the slight time difference, though I didn't measure it to 0.01 seconds :-) What it's doing in that time, I don't know.

The best I can figure is that in FIX 3 it was satisfied with the initial result and presented that. In FIX 4 it wasn't satisfied, but rather than doubling the number of points it first tried some other tweak, which did satisfy it. In FIX 5 the tweak wasn't enough so it does go ahead and double the number of points.

But just what that "tweak" on top of the Romberg extrapolation is, I don't know. And if that tweak only takes half a second, why doesn't it always apply it?

BTW, have you ever looked at all ten digits of a typical error estimate? Where is that other non-zero digit coming from, and why isn't it always in the same place relative to the first one?


Possibly Related Threads...
Thread Author Replies Views Last Post
  Proceedings of the IEEE Article Katie Wasserman 24 2,169 11-29-2013, 07:56 AM
Last Post: Walter B
  XCas / Prime "solve" question Nigel J Dowrick 4 687 11-08-2013, 04:01 AM
Last Post: Nigel J Dowrick
  Using the Prime to solve for eigenvalues Michael de Estrada 28 2,838 10-27-2013, 07:21 AM
Last Post: Tarcisi C
  HP Prime Define: Use with solve, etc. Helge Gabert 0 313 10-23-2013, 06:24 PM
Last Post: Helge Gabert
  HP 34s solve f'(x)=0 Richard Berler 8 747 10-07-2013, 03:03 PM
Last Post: Dieter
  Good puzzle for kids to solve on 35s? snaggs 11 1,077 09-18-2013, 10:40 PM
Last Post: David Hayden
  FYI: HHC 2013 Programming contests are coming Gene Wright 0 318 09-17-2013, 03:12 PM
Last Post: gene wright
  HP-Prime article on Wikipedia CompSystems 4 559 09-08-2013, 06:54 PM
Last Post: Bruce Bergman
  September SOLVE Newsletter Robert Prosperi 4 559 09-07-2013, 01:12 AM
Last Post: Mic
  New article on a new type of neo linear equations Namir 0 387 08-11-2013, 10:27 AM
Last Post: Namir

Forum Jump: