I hope this post does not offend the creator of gsll, as this is not the intent. There are two ways to read this post: One is that it may raise a valid point, and the other that I am missing something already provided, or a mode of using the routines that would make this comment unnecessary.
So, in effect, I think (as a newbie) I am missing something. But what?
I find the naming convention in gsll inconvenient. For example the routine solve-tridiagonal is a wrapper for gsl_linalg_solve_tridiag. To use it to solve a system defined by vectors d e f b, one needs to use a letm block to this effect:
(letm ((d* (vector-double-float d)) (e* (vector-double-float e)) (f* (vector-double-float f)) (x* (vector-double-float (make-array (length d) :element-type 'double :initial-element 0d0)) (b* (vector-double-float b))) (solve-tridiagonal d* e* f* b* x*) (data x*))
Using the routine solve-tridiagonal is still not convenient, because one still needs to convert the argument vectors via the vector-double-float, obscuring the code. If I have multiple calls, I will have to repeat these, or define my own interface function. And here is the problem: *the very nice name solve-tridiagonal is taken*. I could name it solve-tridiag, but sooner or later, I will make a mistake and call the wrong function with the incorrect arguments.
Ideally, once we have gsll, one would just need to do (solve-tridiagonal d e f b). To accomplish that, gsll would need to use the following convention: - the c-code wrapper would be named the same as the c-code, except with the underscores replaced with dashes. Thus we would have (defmfun gsl-linalg-solve-tridiag (diag e f b x) "gsl_linalg_solve_tridiag" (((pointer diag) gsl-vector-c) ((pointer e) gsl-vector-c) ((pointer f) gsl-vector-c) ((pointer b) gsl-vector-c) ((pointer x) gsl-vector-c)) :invalidate (x)... etc)
- Users/contributors can define interface routines that would link to the gsll routine: (defun solve-tridiag (d e f b) (letm ((d* (vector-double-float d)) (e* (vector-double-float e)) (f* (vector-double-float f)) (x* (vector-double-float (make-array (length d) :element-type 'double :initial-element 0d0)) (b* (vector-double-float b))) (solve-tridiagonal d* e* f* b* x*) (data x*))
(I have not tried it, but maybe macros should be better to remove the pass-by-value cost of calling the routine).
So, what am I missing?
Thanks,
Mirko