Re: LAPACK vs. GSL matrix algebra

[Gerard Jungman <jungman at lanl dot gov>]

On Wed, 2004-09-29 at 11:54, Linas Vepstas wrote:
> 
> Hmm. Trying to figure out if I should volunteer to do this.

You would be our hero.

> Google only seems to hit the f2c'ed-up version of lapack.
> Whatever work was done towards that standard seems to be burried
> somewhere.

The interfaces expressed in the ATLAS clapack.h might be
(part of) what I am thinking of. Clint Whaley must have some
notion of what to do for a complete wrapper, although
the ATLAS clapack.h interfaces have only some of the
functions. I vaguely remember some discussion about
standardization of C interfaces years ago. The C-BLAS
interface draft came out around 1998, and there may
have been some discussion of lapack interfaces at that
time.

> I could try to hand-build a more acceptable binding, but would
> like to have active review for that.

Right. The logical thing would be to create a generic C
interface and then target GSL to that. That means you
would effectively become the C-LAPACK standardizer.
We would need to get some input on that, probably at
least from Whaley if nobody else.

I guess the main issue is handling the workspace stuff,
although I haven't thought about it. It needs to be simple
but flexible. Brian thought about this about 5 years ago,
so he should say something...

> Bing.  So I guess that you are thinking of just literally 
> copying c-lapack into the gsl cvs tree, and using that as 
> the starting point?

It is probably necessary to have a "native" GSL implementation,
in the same way that we have a native cblas, just so people
have something that works out of the box. I originally wanted
to just dump a fortran blas implementation into a subdirectory
of GSL and use the cblas reference implementation to provide
a cblas; nasty and inefficient, but hardly any worse than our
hand-rolled cblas. Brian argued that we should roll our own,
I think based on the grodiness of having to carry around
and build fortran sources as an integral part of GSL.

Copying clapack into the source tree is one way to get the
functions to build out of the box. But I'm open to any ideas.
Also, it doesn't solve the general problem of creating
reasonable C interfaces; we would probably have to do
that anyway.

I wish that ATLAS had full lapack functionality. Then at
least we would have that to look at. I don't know what
his plans are in this regard.

> Are there any license problems? lapack seems to be either 
> bsd licnese, or public domain; I can't tell which; gsl is gpl'ed.

lapack itself is certainly free enough. The rpm that I have was
built by Redhat, from the netlib sources, and the "License" entry
says "Freely distributable". I couldn't find any markings on
the f2c'ed clapack sources, but I'm guessing that they can't
be more restrictive than lapack and f2c themselves.

> The "professional" way to solve this is to do it at run-time, 
> by calling through a table of function pointers.  By default,
> the table points to gsl-cblas, but either the whole teable or
> individual routines in it can be replaced by the user, as desired.

I guess that makes sense. It makes me a little nervous.

> Or maybe the function pointer table would be initialized 
> with atlas routines I guess ... 

In principle the libraries all use the standard cblas_XXX names.
How do you get all these linked together so that they are
available together at run time? Is this a dumb question?
This is the main reason that I thought the link-time
solution was natural.

> Do you want a sketch of what the function-pointer table might 
> look like?

Sure. It needs to be thread-safe, but that's probably not hard.
How do you solve the shared name problem?

-- 
Gerard Jungman <jungman@lanl.gov>
Los Alamos National Laboratory