Estimating Accumulated Rounding Errors


Hello all,

I am working on yet another Pi program. In the past to eliminate accumulated rounding errors I just added extra guard digits. However I'd like to be able to precompute exactly how many I need. I can with great certainty compute the number and type of operations, and without much effort I can produce a table at regular intervals of the number of required guard digits. But, I'd like to explore a pure mathematical approach.

I have a few ideas, but I'd also like your thoughts as well.


Edited: 15 Aug 2012, 12:45 a.m.


Start with Message #24 & enjoy reading lots :)

- Pauli


Pauli, thanks. I think I have about 1/3-1/2 of those books.


You might find error discussion in the appendix of the HP-15C Advanced Functions Handbook of some use. (And certainly of more use than my thoughts on the topic!)


don't estimate the errors - know them! or rather their upper bounds.

i suggest the field of validated computing, which i have been looking at for
some time. This is the use of interval arithmetic to bound your computations.

you have to define a new number class with the usual operators, which holds an
interval [Inf, Sup]. Where your true real value is somewhere Inf <= x <=
Sup. The operators maintain the error window whenever you perform
arithmetic. If you have, say a Taylor series of arctan, you have to also
expand the window by the max error of your truncation and so on.

some years ago i successfully built a calculator around this idea.
(windows exe)


exact -acc 2000
increasing precision to 2012 =

which automatically increases precision until the error bound is met.

Back then i was planning to port this onto physical calculators like the
hp50g, but i decided i didn't need silly accuracy but only enough.

However, i've always been tantalised by the idea of applying the same approach
to solving once and for all the calculator problems of solve and
integrate. The problem of bounding the error with classical algorithms fails
because they are all based on sampling.

A rigorous numerical approach works for numerical quadrature by dynamically
creating a taylor expansion of your function (numerically) and essentially
integrating term by term, whilst keeping check on the overall error
bounds. it's like mixing numeric with symbolic - kind of semi-symbolic where
you only have to manipulate a series.

Anyhow, progress has been made in this area since i last looked. just
recently, i've been reading "Validated Numerics by Warwick Tucker." This is
quite a good introduction to the topic.

Edited: 15 Aug 2012, 8:28 p.m.


weird, everyone else' posts wrap except mine.


The editor is too polite to break you pi :)


If you enclose something in [pre] .... [pre/] then the text will remain untouched, no automatic line break is done. You need to do the line feeds manually, as I've done below. As I didn't want to edit your post, you may replace your one long line with a copy of my reformatted lines if you wish, in order to prevent us from having always long lines throughout the entire thread.


exact -acc 2000
increasing precision to 2012 =



i did [pre] .. [/pre]

always making this mistake...


I too have thought about interval arithmetic in a calculator :-)

I was specifically thinking of the mpfi library which is based on the guarantee correctly rounded MPFR library. They are not decimal arithmetic however.

- Pauli


Hi Paul,

Regarding interval arithmetic, i've been experimenting with C-XSC. but this is mainly for researching the subject.

For calculators, i'm thinking of extending my own number library. There's a neat way to encode an interval without much overhead:

intervals have a "mid" (m) and a "radius" (r) so [a,b] is the same as [m-r,m+r]. so you store the mid value in your float format then encode the radius as the last digit of the mantissa. Naturally, you have to know in your implementation not to use this last digit as part of the mid value. but the clever idea is you don't really need much accuracy in the radius, which now has a limited mantissa (eg 1 digit) and its magnitude correctly lies just below that of the mid (since otherwise you're wasting mid mantissa that overlaps the error anyhow).

Actually, just 1 digit of error term is a bit poor, the consequence is wider error bounds and unnecessary extra calculation. However, you may recall that my decimal library works with "digits" of 4 decimal digits, ie base 10,000. I would dedicate one of these such digits and benefit from a 4 digit (decimal) error mantissa. I'm thinking this could be an unexpected bonus.

Regarding multiple precision, there are two possibilities here; either make my library extended precision, which will cost me another "digit" for the length and therefore take up more space and time. or, stay with my 25 digits of working precision and just allow the error to swallow what it takes and hope to deliver enough validated payload for people's actual need (eg 10 say). I think if it could deliver 10 validated digits of numerical quadrature, that would be a big result.

i'm also bearing in mind that, although calculators are much more powerful now than ever, it can't be completely crazy.

-- hugh.


A couple of thoughts/questions.

Wouldn't your proposed scheme for m, r using the lowest digit block limit the width of the intervals? I'd be rather surprised if intervals didn't get larger than any finite limit on the width. i.e. high accuracy isn't required in the width, but large range is.

It also seems like there might be rounding issues determining the bounds of the interval -- perhaps depending on how it were implemented. This could increase the effective width unnecessarily.

Still interval arithmetic is an interesting idea and one I'd like to have included in the 34S :-)



You might be right for certain cases. This is something i have to look into. For example, sometimes there's need for things like [x,Inf] or [-Inf,x] to represent upper or lower possible values. What i don't know is whether, such cases come out in normal calculation AND lead on to useful results. If it means your result has gone wide or these are simply means to encode an exception or divide by zero, then perhaps we don't really need them.

my working idea is that the range, when it encodes uncertainty, is within the magnitude of the value. Otherwise your value is useless. so when x = 0.5 with uncertainty +/- 10, the 0.5 isnt really telling us anything.

if this idea is workable, i can sneak the uncertainty into the number representation with fairly little overhead. I admit it's not going to encode ranges greater than the number itself. is that a problem; i don't know yet.

on the subject of rounding errors, the current fashion for implementation is to use a pair of IEEE doubles for your interval [a,b], and when `a' is calculated you set the floating point mode to round DOWN and when `b' is calculated, you set the float mode to round UP.

you have to mess about with setting the float mode (which is usually not very portable), but it means you can use your built-in hardware floating point, which is usually present. This is quite a neat idea and you're sure to maintain a hard bound.

incidentally, to take an example; what they do when working on binary floats (eg almost always) with numbers like 0.1 is to input "[0.1,0.1]" so they get the greatest representable float <= 0.1 and the least representable float >= 0.1. in other words even some "constants" wind up being a non-zero interval.

-- hugh.


Hugh, thanks. I think this is what I was looking for.

Possibly Related Threads…
Thread Author Replies Views Last Post
  HP Prime program: rounding to a fraction Patrice 3 1,709 10-31-2013, 06:16 AM
Last Post: Joe Horn
  HP Prime: Rounding error in determinant Stephan Matthys 3 1,647 10-25-2013, 09:29 PM
Last Post: Walter B
  HP Prime editor errors Alberto Candel 0 972 10-22-2013, 01:34 AM
Last Post: Alberto Candel
  Rounding of 10^x Olivier De Smet 8 2,557 08-28-2013, 06:33 AM
Last Post: Dieter
  Paper on the errors in HP Calculators vs some older computer processors Les Koller 3 1,773 08-19-2013, 12:37 PM
Last Post: Mike Morrow
  Estimating accuracy in finite precision computations mpi 17 4,350 02-22-2013, 09:44 AM
Last Post: mpi
  HP-28S Syntax Errors William N Strew 7 2,562 07-15-2012, 05:09 PM
Last Post: Luiz C. Vieira (Brazil)
  Some errors about NOMAS, (c) of 71B ROM images Christoph Giesselink 4 1,791 06-13-2012, 04:10 AM
Last Post: Valentin Albillo
  Errors reading with a HP 82153A Barcode Optical Wand aurelio 22 5,886 04-04-2012, 05:19 PM
Last Post: Les Wright
  WP34s: number display and rounding Dieter 4 1,459 10-14-2011, 05:11 PM
Last Post: Marcus von Cube, Germany

Forum Jump: