Re: [Gsll-devel] small modification to blas2.lisp

On Wed, May 14, 2008 at 2:23 PM, Mirko Vukovic mirko.vukovic@gmail.com wrote:

Hi,

I was playing with gsll and blas. I modified the cblas-transpose definition in blas2.lisp (see top of file) to capture a few more values:

Here is the new version:

(cffi:defcenum cblas-transpose "CBLAS_TRANSPOSE from /usr/include/gsl/gsl_cblas.h." (:notrans 111) (:trans 112) (:conjtrans 113))

I fail to see the difference with what is in there now: (cffi:defcenum cblas-transpose "CBLAS_TRANSPOSE from /usr/include/gsl/gsl_cblas.h." (:notrans 111) :trans :conjtrans) See http://common-lisp.net/project/cffi/manual/html_node/defcenum.html#defcenum "If value is omitted its value will either be 0, if it's the first entry, or it it will continue the progression from the last specified value." so :trans gets 112, :conjtrans 113 by default.

I am posting it here because I was not able to download with svn so I don't have versioning for gsll.

Also, here is a code snippet that I used to test for tridiagonal matrix solving that also uses gsll-s blas2 routine gemv to double check that things work.

(in-package :gsll) ;; solve a tridiagonal system and then verify result by matrix multiplication ;; uses solve-tridiagonal and gemv (let* ((dim 6) (dim- (- dim 1))) (let ((d (make-array dim :element-type 'double-float :initial-element -2.d0)) (e (make-array dim- :element-type 'double-float :initial-element 0.2d0)) (f (make-array dim- :element-type 'double-float :initial-element 1.d0)) (x (make-array dim :element-type 'double-float :initial-element 0.d0)) (b (make-array dim :element-type 'double-float :initial-element 0.5d0)) (mat (make-array (list dim dim) :element-type 'double-float :initial-element 0.d0)))

;; load matrix
(loop for i from 0 below dim-
   do (progn
        (setf (aref mat i i) (aref d i))
        (setf (aref mat (+ i 1) i) (aref f i))
        (setf (aref mat i (+ i 1)) (aref e i))))
(setf (aref mat dim- dim-) (aref d dim-))

;; gsl operations
(letm ((d* (vector-double-float d))
       (e* (vector-double-float e))
       (f* (vector-double-float f))
       (x* (vector-double-float x))
       (b* (vector-double-float b))
       (v2* (vector-double-float
                   (make-array dim :element-type 'double-float
                                   :initial-element 0.d0)))
       (mat* (matrix-double-float mat)))
  (solve-tridiagonal d* e* f* b* x*)
  ;; gemv returns the result in the last argument
  (gemv :notrans 1.d0 mat* x* 0.d0 v2*)
  (data v2*))))

I am sure this code betrays my lack of lisp experience. So feel free to clean it up if you want it included in the verification part.

Mirko

I might make a few simplifications; #'1- to do decrementing and I don't think there's a need to make arrays in CL and then remake them on the foreign side. It is possibly a candidate for example/ regression tests.

Liam