From: rcarter@gaudi.nas.nasa.gov (Russell L. Carter)
Subject: Examples anyone? was Re: 128 bit floating point
Followup-To: sci.math.num-analysis
Organization: NAS, NASA Ames Research Center, Moffett Field, CA
Date: Fri, 17 Apr 92 16:12:44 GMT

In article <36358@ogicse.ogi.edu> borasky@ogicse.ogi.edu (M. Edward Borasky) writes:

[quoting my original post]

>>I'd love to hear about particular examples of algorithms that
>>really need 128 bit fp and don't have more accurate alternatives.

and then he comments:

>The classical example is the solution of ill-conditioned linear equations.
>The standard algorithm (LINPACK, LAPACK) uses Gaussian elimination (LU
>decomposition) with partial pivoting to factor the matrix.  This takes 
>approximately 2/3*N**3 flops (adds AND multiplies) and its error 
>properties are well known ("Matrix Computations" by Golub and Van Loan).
>For a given type of matrix, eventually as N gets large enough the 
>condition number will get bad enough that LU decomposition in 64 bits no
>longer give you a usable answer, and you must go to 128 bits. 

Your response is theoretical, and I am fully aware of the theory.
My question still stands, and I would now *really* like to hear
of some examples that differ from the contrived ones of Kahan.
In particular, I would love to find examples arising in the
simulation of physical systems, and I would writhe in ecstasy
if that example happened to require the solution of a really
large scale pde problem.

After recently hearing the story about how native 128 bit
floating point infected a certain forthcoming high performance 
microprocessor design, I am becoming convinced that religion more 
than anything else is driving this issue.

>
>But there is another algorithm, based on QR decomposition rather than
>LU, which doesn't have the numerical problems that Gaussian elimination
>exhibits (also in Golub and Van Loan).  This algorithm requires twice 
>the flops as LU.  Unfortunately, the number of multiplies and adds in 
>the inner loop is not the same; most supercomputers do one of each per 
>cycle, so you won't get peak speed on a QR decomposition.  On the LU, 
>whether you use dot product or "GAXPY" you have one add and one multiply

Gee, if the matrix is full, it seems like using SGEQRF from LAPACK
(90% efficient on 1 cpu, >70% efficient on 8 cpus of a Y-MP) would
be appropriate.  If the matrix is sparse AND has poor conditioning
AND it arose from trying to solve a real physical problem, and
most importantly, there's no other practical way to solve the problem,
then don't be shy, tell me about it!  I'd really like to find
one of these.

BTW, it is interesting that the same folks that brought us LAPACK
also proved that in general, sparse system solves have exceptionally
nice error behavior.  The error refinement work is very impressive,
and I think it will probably be the most enduring result of the project.

>in each cycle; the computation runs at near-peak speed in today's 
>machines.  In fact, special-purpose machines (FPS-164/MAX) exist that 
>do the LU decomposition.

Evidently, there must not be much of a market.

>
>My question, then, is, given all this, why don't those who have large
>poorly-conditioned linear equations to solve build or sponsor special-
>purpose machines based on the QR decomposition?  It's been proposed,
>but never, to my knowledge, been done.

Why build something special purpose when the general purpose system
runs at 90% efficiency?  You can have your cake and emacs, too.

Russell Carter
Parallel Systems Development
Numerical Aerodynamic Simulation Facility

rcarter@nas.nasa.gov
Newsgroups: comp.arch