GSL v2.0 discussion

[Patrick Alken <patrick dot alken at Colorado dot EDU>]

Hello all,

  There have been a lot of new features added to GSL since the last 
release 1.16, some of which have modified underlying data structures 
(breaking binary compatibility) and even some API changes. For this 
reason I don't want to make another 1.x release, but I know at the same 
time many people were hoping that 2.0 would interface with advanced 
linear algebra libraries (like libflame or lapack).

My current todo list before the next release is:
1) import interp2d code (mostly done just needs some further testing)
2) Hermite interpolation (shouldn't take long once I find some time)

Adding libflame or lapack to the todo list would delay the next release 
significantly. But I want to open this discussion again (it has been 
discussed several times previously on this list).

One method of doing this is to make 3 new libraries:
1) libgsllinalg - contains all the current GSL implementations of the 
linear algebra routines
2) libgslflame - contains wrapper routines with the same interfaces that 
gsl_linalg has now, but calls the relevant libflame routine
3) libgsllapack - contains wrapper routines with the same interfaces of 
gsl_linalg, but calls the relevant lapack routine

This way the user can decide which linear algebra library they want at 
link time, while maintaining a common API for everything. If the user 
wants to use the current GSL routines, they would link with -lgsl 
-lgslcblas -lgsllinalg

If they want libflame, they would link with -lgsl -lgslcblas -lgslflame 
-lflame

If they want lapack, they would link with -lgsl -lgslcblas -lgsllapack 
-llapack

A second method of doing this is to decide at compile time which linear 
algebra library to use (ie: autoconf could detect if lapack or flame is 
installed and then compile the library accordingly making the correct 
calls to the lapack/flame routines). This has the advantage of keeping 
the -lgsl -lgslcblas link libraries the same, the user would just need 
to add -llapack or -lflame.

*Potential issues*

1. Workspace sizes

Many of the flame/lapack routines require workspaces of different sizes 
than each other and the corresponding GSL routines. The easiest way to 
handle this would be to dynamically allocate workspace inside the 
wrapper functions of the gslflame and gsllapack libraries, which is not 
too appealing.

Another way to handle it is to have gsl_linalg_alloc() routines for all 
functions which would figure out the right workspace size for whatever 
library you are using. This would require a lot of code changes for 
current users, also not too appealing.

2. Errors

LAPACK tends to handle errors pretty well, and returns error codes from 
each function which could be interpreted by the GSL wrapper and passed 
on to the user.

libFLAME however seems to lack an error handling mechanism and calls 
abort() when it encounters an error (I have not looked deeply into the 
source but this is based on previous discussions on this list). Calling 
abort() is a very poor way to handle this and GSL users could be 
frustrated by this. I don't want to get into the business of 
forking/rewriting the flame error handling stuff so I don't see any easy 
solution to this.

3. Global state

This is another advantage of LAPACK over libflame - libflame requires a 
global initialization call to initialize a small number of global 
variables. This means that a) each wrapper function for libflame would 
need to check if that memory is initialized and if not call the function 
or b) GSL needs to provide a global init/free call. I don't like either 
of these options. (a) is bad because there is no simple way for the user 
to free that memory from their program. (b) is bad because global states 
are just ugly and bad design imo.

4. Fortran vs C

One main disadvantage of LAPACK is the ordering of memory in fortran, 
which would require transposing matrices prior to the lapack call, and 
then transposing back afterward. This isn't such a big problem in terms 
of efficiency but its still a headache.

*Summary*

I would like feedback from everyone on the following points:

1. Should we try to add lapack/flame interfaces for the 2.0 release or 
wait until 3.0? I personally probably won't have a lot of time to work 
on this for several months.

2. Is it better to select gsllinalg/lapack/flame at compile time or link 
time?

3. Whats the best way to handle the memory workspace requirements (add 
_alloc functions to all gsl_linalg routines or dynamic allocation)?

4. What should we do about error handling in libflame? Just accept the 
abort() behavior?

5. Is there a strong preference for doing wrappers for both lapack and 
flame? Should we only interface to lapack, due to the difficulties with 
flame (global state, abort() error handling)? Should we only interface 
to flame due to its more modern design?

Patrick