This site is updated infrequently. For up-to-date information, please visit the new OCaml website at ocaml.org.

RE: [Caml-list] Efficient and canonical set representation?
[ Home ] [ Index: by date | by threads ]
[ Search: ]

[ Message by date: previous | next ] [ Message in thread: previous | next ] [ Thread: previous | next ]
 Date: 2003-11-10 (20:36) From: Eric Dahlman Subject: Re: [Caml-list] Rounding mode
```Brian Hurt wrote:

>On Mon, 10 Nov 2003, Eric Dahlman wrote:
>
>
>
>>Christophe Raffalli wrote:
>>
>>
>>
>>>I saw a previous discussion about rounding mode for OCaml (in 2000).
>>>
>>>As anyone implemented the functions to change the rounding mode for
>>>floating point from OCaml (possibly as a patch to add primitives and
>>>save a C function call) ?
>>>
>>>This is necessary to implement interval arithmetic ...
>>>
>>>
>>Somewhat off topic but why is this necessary from a numerical math type
>>of perspective.  I am honestly curious as I don't see how this would
>>interact with the calculation in a meaningful way.
>>
>>
>>
>
>Simplistic discussion of why changing the rounding modes is usefull for
>numeric programming.  If you know numerical analysis, stop reading now
>:-).
>
>
I read it anyway ;-)

[Snip a description of how floating point works.]

>Consider, for example, multiplying two numbers- x and y.  But you don't
>have x and y exactly, you have (x + e) and (y + e) for some error amount
>e.  You multiply them, and you get (x + e)*(y + e) = x*y + x*e + y*e +
>e*e.  Now, e*e is vanishingly small, we can ignore it.  And, assuming x
>and y are both about 1.0, our answer is x*y + ~2*e- in other words, we've
>about doubled the amount of error.  Even starting at 1e-15, doubling it
>every operation and you eat up all your precision surprisingly fast.
>
>
>How can we tell how much error we have?  Well, you can either do a whole
>bunch of graduate-level mathematics, or you can run the same algorithm
>with different rounding modes.  If you always round towards positive
>infinity, you'll always be adding an error term- i.e. it'll be (x + e) and
>(y + e) with a result of x*y + 2e.  If you always round to negative
>infinity, you'll be subtracting the error term- (x - e) and (y - e) with a
>result x*Y - 2e.  Comparing the results gives you both a good idea of
>where the answer really is, and how big your errors are.
>
>
Now I don't know that it is this simple.  I understand that by always
rounding in one direction you are in effect  calculating a bound on the
error in that direction but that will only work in the most simplistic
of cases.  For instance if half of your values had their sign changed
before being summed the sign of half  the associated error terms would
also change.  The result would not give you any information about the
effective lower/upper bound you thought you were calculating.  It seems
to me that this approach is only applicable to problems which would be
trivially solvable with graduate level mathematics and of absolutely no
value in any interesting cases.

>Note, all is not necessarily lost.  If instead we have (x + e) and (y - e)
>a little bit of quick algebra shows that the result is actually *more*
>accurate than either of the original inputs.  The errors tend to cancel
>instead of accumulate.  Numerical analysts call algorithms whose errors
>tend to cancel "stable", while algorithsm whose errors tend to accumulate
>are "unstable".  A good textbook on numerical methods will generally only
>discuss stable algorithms- use 'em.
>
>
This is incorrect because you are thinking about real numbers while
working with floating point numbers, this is especially dangerous when
you are operating near or below the epsilon for the particular floating
point representation like we are now.  If you were to look at the set of
numbers which are precisely represented by a given floating point
representation you would see that it is concentrated near 0 and that it
tends to thin out as you move away from there.   The distance  between
any  two adjacent floating point numbers is highly variable and not
uniform or "nice" by almost any measure.  I really am having a hard time
figuring out what you can say about the tiny e terms in this case as we
move between regions with differing "floating point resolutions".

Just to illustrate this the different "resolutions" the smallest
positive double float which can be represented is  2.225 x 10^(-308)
while the epsilon for doubles in general is 1.11 x 10^(-16) and can be
thought of at the coarsest resolution.  That is a big difference and is
because floats are packed in tightly around 0.

As for a numerically stable algorithm being able to produce a result
which is more accurate than the original inputs that is total hogwash.
You cannot ever have more accuracy than was present in your original
input, ever. Nihil ex nihilo.

>The effect can also be sort of acheived by using a random rounding method.
>For each operation, the CPU picks one of it's rounding methods at random
>to use (well, pseudo-random anyways).  Now, while this will in general
>increase the accuracy of floating point programs, it also means that every
>time you run a numerical program, even over the same data, you will get a
>slightly different answer.  Newbies wondering why their program keeps
>changing it's answer will probably rival newbies wondering why their stack
>overflows pretty quickly.  Which is why I'm hesitant to advocate doing
>this to all programs.  I still would like to be able to select the
>rounding mode- either at compile time or at run time, which ever is
>easier.  Actually, run time is probably easier to implement and more
>usefull (I can switch rounding modes for different parts of the program).
>
>
As an oldby I would be disconcerted by this as well.  All you are doing
is trading one sort of  numerical noise for another which you guess
should be better behaved but in the absence of a neat proof to the
contrary I don't believe is.  There are just *SO* many sources of huger
error terms in numerical calculations that the tiny rounding errors are
really lost in the noise.  I even tried to come up with a case where it
would matter but that got so farcical that I gave up on it.

What I was really wondering about when I first asked this question was
if the other person had exactly that kind of example where they were
using some clever trick to usefully misapply floating point calculations
in their problem domain.

-Eric

-------------------
To unsubscribe, mail caml-list-request@inria.fr Archives: http://caml.inria.fr
Bug reports: http://caml.inria.fr/bin/caml-bugs FAQ: http://caml.inria.fr/FAQ/
Beginner's list: http://groups.yahoo.com/group/ocaml_beginners

```