blas.lisp (3680B)
1 (defpackage :blas 2 (:use :common-lisp :sb-alien :sb-ext :sb-sys) 3 (:export :blas-float 4 :+blas-float-size+ 5 :+blas-0+ 6 :+blas-1+ 7 :+blas-i+ 8 :*blas-alien-0* 9 :*blas-alien-1* 10 :make-blas-array 11 :make-integer-array 12 :make-double-array 13 :blas-array-alien 14 :integer-array-alien 15 :double-array-alien 16 :dlamch 17 :zaxpy 18 :zgemm 19 :zgesvd 20 :zgesdd 21 :zgesvdx 22 :zheevr)) 23 24 (in-package :blas) 25 26 (deftype blas-float () '(complex double-float)) 27 28 (defconstant +blas-float-size+ 16) 29 30 (defconstant +blas-0+ (coerce 0 'blas-float)) 31 (defconstant +blas-1+ (coerce 1 'blas-float)) 32 (defconstant +blas-i+ (coerce (complex 0 1) 'blas-float)) 33 34 (defparameter *blas-alien-0* 35 (make-array 2 :element-type 'double-float :initial-contents '(0d0 0d0))) 36 37 (defparameter *blas-alien-1* 38 (make-array 2 :element-type 'double-float :initial-contents '(1d0 0d0))) 39 40 (defmacro make-blas-array (dimensions &rest args) 41 `(make-array ,dimensions ,@args :element-type 'blas-float)) 42 43 (defmacro make-integer-array (dimensions &rest args) 44 `(make-array ,dimensions ,@args :element-type '(signed-byte 32))) 45 46 (defmacro make-double-array (dimensions &rest args) 47 `(make-array ,dimensions ,@args :element-type 'double-float)) 48 49 (defun blas-array-alien (array) 50 (sap-alien (vector-sap (array-storage-vector array)) (* double))) 51 52 (defun integer-array-alien (array) 53 (sap-alien (vector-sap (array-storage-vector array)) (* int))) 54 55 (defun double-array-alien (array) 56 (sap-alien (vector-sap (array-storage-vector array)) (* double))) 57 58 (declaim (inline dlamch zaxpy zgemm zgesvd zgesdd zgesvdx zheevr)) 59 60 (define-alien-routine ("dlamch_" dlamch) double 61 (cmach c-string)) 62 63 (define-alien-routine ("zaxpy_" zaxpy) void 64 (n int :copy) 65 (za (* double)) 66 (zx (* double)) 67 (incx int :copy) 68 (zy (* double)) 69 (incy int :copy)) 70 71 (define-alien-routine ("zgemm_" zgemm) void 72 (transa c-string) 73 (transb c-string) 74 (m int :copy) 75 (n int :copy) 76 (k int :copy) 77 (alpha (* double)) 78 (a (* double)) 79 (lda int :copy) 80 (b (* double)) 81 (ldb int :copy) 82 (beta (* double)) 83 (c (* double)) 84 (ldc int :copy)) 85 86 (define-alien-routine ("zgesvd_" zgesvd) void 87 (jobu c-string) 88 (jobvt c-string) 89 (m int :copy) 90 (n int :copy) 91 (a (* double)) 92 (lda int :copy) 93 (s (* double)) 94 (u (* double)) 95 (ldu int :copy) 96 (vt (* double)) 97 (ldvt int :copy) 98 (work (* double)) 99 (lwork int :copy) 100 (rwork (* double)) 101 (info int :out)) 102 103 (define-alien-routine ("zgesdd_" zgesdd) void 104 (jobz c-string) 105 (m int :copy) 106 (n int :copy) 107 (a (* double)) 108 (lda int :copy) 109 (s (* double)) 110 (u (* double)) 111 (ldu int :copy) 112 (vt (* double)) 113 (ldvt int :copy) 114 (work (* double)) 115 (lwork int :copy) 116 (rwork (* double)) 117 (iwork (* int)) 118 (info int :out)) 119 120 (define-alien-routine ("zgesvdx_" zgesvdx) void 121 (jobu c-string) 122 (jobvt c-string) 123 (range c-string) 124 (m int :copy) 125 (n int :copy) 126 (a (* double)) 127 (lda int :copy) 128 (vl double :copy) 129 (vu double :copy) 130 (il int :copy) 131 (iu int :copy) 132 (ns int :out) 133 (s (* double)) 134 (u (* double)) 135 (ldu int :copy) 136 (vt (* double)) 137 (ldvt int :copy) 138 (work (* double)) 139 (lwork int :copy) 140 (rwork (* double)) 141 (iwork (* int)) 142 (info int :out)) 143 144 (define-alien-routine ("zheevr_" zheevr) void 145 (jobz c-string) 146 (range c-string) 147 (uplo c-string) 148 (n int :copy) 149 (a (* double)) 150 (lda int :copy) 151 (vl double :copy) 152 (vu double :copy) 153 (il int :copy) 154 (iu int :copy) 155 (abstol double :copy) 156 (m int :in-out) 157 (w (* double)) 158 (z (* double)) 159 (ldz int :copy) 160 (isuppz (* int)) 161 (work (* double)) 162 (lwork int :copy) 163 (rwork (* double)) 164 (lrwork int :copy) 165 (iwork (* int)) 166 (liwork int :copy) 167 (info int :out))