Re: complexity of numerical software



David N. Williams wrote:
This is a repost of a message which appeared at my news server,
but hasn't showed up in Google groups for several hours.  My
apologies if it's actually duplicated.

Jean-Claude Arbaut wrote:

Jaap Spies wrote:


Jean-Claude Arbaut wrote:


Jaap Spies wrote:


David N. Williams wrote in a reaction to carlos@xxxxxxxxxxxx:


I'm really puzzled by your classification of numerical software.
I would have put it in the truly challenging category!?

-- David


Numerical software? What kind of numerical software you are talking
about?


Maybe it would be intersting to know your definition of "numerical
software"...


Yes. We have the same question!


Okay. I thought we were talking about code like SLATEC, NSWC, CMLIB,
etc. They have (almost) nothing to do with algorithms to find subsets,
which I would classify with discrete or combinatorics algorithms.


SLATEC, etc., are indeed examples of the kind of thing I had in
mind.


It's easy to write poor numerical algorithms (with my definition),
the usual example is "Numerical Recipes", but I find much more
difficult to write *good* (read /accrurate and fast/) ones. Just
try to write a correctly rounded sine (or even a square root !).


Jaap Spies also wrote:

[...]

Testing and debugging of numeric algorithms is easy. Just check the
results.


Actually, I would include testing numerical analysis software in
the truly challenging category.  A classic example is William
Cody's celefunt suite of accuracy tests for complex elementary
functions.  It makes clever use of mathematical identities, but
while useful it is still limited, as Cody points out in his ACM
article on Algorithm 714.  Surely it is nontrivial in general to
analyze where inaccuracies for a numerical program are likely to
occur, to even begin testing.

Jean-Claude mentioned the square root.  Here's a quote (up to
math symbols) from one of those truly wonderful old IBM manuals
for VS FORTRAN:

  Effect of Argument Error

  epsilon ~ delta/2

  Accuracy

  The accuracy of the algorithm is perfect; for all (15*2^27)
  nonnegative normalized short floating-point numbers x, SQRT
  returns the correctly rounded value...

There's clearly nontrivial error analysis underlying IBM's
implementation of SQRT,

I saw that in MacOSX implementation of function sqrt, which I suspect comes from IBM: two Newton iteration, plus a correction on the last bits. Everybody can write a Newton iteration, but the remaining part is what makes good algorithms...

which you have to worry about if you
don't have access to a trusted implementation like theirs.  Then
the known relationship between the input error, delta, and the
output error, epsilon, maybe gives you some hope of
understanding the results of accuracy tests on code that uses
SQRT.

And the subtleties of floating-point arithmetic feed in as well.

-- David

Very difficult to write, indeed ;-) .



Relevant Pages

  • Current standings: Structures
    ... David N. Williams ... SFIELD: 2 cells ... The semantics of DFIELD: ...
    (comp.lang.forth)
  • Re: matrix convolution and fourier transform
    ... Dave Robinson wrote: ... > David Epstein wrote: ... > for explaining what its about, and not trying to impress you by ... > uses the fast frequency domain algorithms, ...
    (comp.soft-sys.matlab)
  • Re: Memoir Class: Title and Chapter/Section on Same Page
    ... > David N. Williams wrote: ... LaTeX FAQ: http://www.tex.ac.uk/faq ...
    (comp.text.tex)
  • Re: Doctors Advice?
    ... Would you like me to offer some posting references and search criteria ... Williams; watch my lips you silly old bastard: ... What perplexes me is that he has created a troll faq where there are 4 ... David, don't you realize that he simply stumble across ...
    (misc.fitness.weights)
  • Re: Receiving error message
    ... On Mon, 2004-08-23 at 13:16, Williams, David wrote: ... I think in subsequent updates this was fixed, at least on one box I ... But I think a recent update renamed the XF86config file and it now uses ...
    (Fedora)