Artur:
Most of the BLAS 1 are hiding in Fortran 90:
*DOT is DOT_PRODUCT,
I*AMAX(A) is MAXLOC(ABS(A)) -- I wanted MAXABSLOC because I don't
trust compilers to optimize away the temp resulting from ABS(A).
*AXPY can be written directly as A*X + Y where A is scalar and X and Y
are conformable.
*ASUM(A) can be written SUM(ABS(A)) (same remarks about temp)
*SCAL can be written directly as A*X
*COPY can be written directly as X = Y
Missing are *SWAP, *ROT*, *NRM2.
It's a "quality of implementation" issue whether DOT_PRODUCT on single-
precision reals is SDOT or SDSDOT, whether D = DOT_PRODUCT(R1, R2)
uses DSDOT or SDOT (D is double precision, R1 and R2 are single), and
whether the compiler optimizes away the temp in DOT_PRODUCT(CONJG(X),Y)
(which you must write to get CDOTU, not CDOTC).
One _could_ insure getting DSDOT by writing SUM(DPROD(R1, R2)) but the
same remarks about the apparent temp apply.
One _could_ compute *NRM2(X) as SQRT(DOT_PRODUCT(X,X)) but *NRM2 is
probably more efficient, doesn't underflow, and doesn't overflow unless
the final result does.
Does MATMUL cover several of the BLAS 2 routines, e.g. GEMM?
By the way, the original papers on BLAS 1 were written by my colleagues
at JPL. Charles Lawson has retired. Richard Hanson left JPL decades ago,
and is now the only remaining mathematician at IMSL (actually, only half-
time -- how serious are IMSL about the "M" part of their name?). Dave
Kincaid was a post-doc, here for just one year. Fred Krogh is all that's
left of the internal computational mathematics support group, and he's
planning to retire if he gets any more administrative headaches. I guess
when Fred retires, we'll do our orbit determination in Excel (it would be
better than Numerical Recipes, anyway).
Best regards,
Van Snyder
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
|