source: project/release/4/blas/trunk/blas.scm @ 26338

Last change on this file since 26338 was 26338, checked in by Ivan Raikov, 8 years ago

blas: using bind instead of easyffi

File size: 79.3 KB
Line 
1;; blas.scm
2
3(module blas
4
5        (RowMajor
6         ColMajor
7         NoTrans
8         Trans
9         ConjTrans
10         Left
11         Right
12         Upper
13         Lower
14         Unit
15         NonUnit
16
17         sicopy
18         dicopy
19         cicopy
20         zicopy
21
22         scopy
23         dcopy
24         ccopy
25         zcopy
26
27         unsafe-sgemm!
28         unsafe-dgemm!
29         unsafe-cgemm!
30         unsafe-zgemm!
31         sgemm!
32         dgemm!
33         cgemm!
34         zgemm!
35         sgemm
36         dgemm
37         cgemm
38         zgemm
39
40         unsafe-ssymm!
41         unsafe-dsymm!
42         unsafe-csymm!
43         unsafe-zsymm!
44         ssymm!
45         dsymm!
46         csymm!
47         zsymm!
48         ssymm
49         dsymm
50         csymm
51         zsymm
52
53         unsafe-chemm!
54         unsafe-zhemm!
55         chemm!
56         zhemm!
57         chemm
58         zhemm
59
60         unsafe-ssyrk!
61         unsafe-dsyrk!
62         unsafe-csyrk!
63         unsafe-zsyrk!
64         ssyrk!
65         dsyrk!
66         csyrk!
67         zsyrk!
68         ssyrk
69         dsyrk
70         csyrk
71         zsyrk
72
73         unsafe-cherk!
74         unsafe-zherk!
75         cherk!
76         zherk!
77         cherk
78         zherk
79
80         unsafe-ssyr2k!
81         unsafe-dsyr2k!
82         unsafe-csyr2k!
83         unsafe-zsyr2k!
84         ssyr2k!
85         dsyr2k!
86         csyr2k!
87         zsyr2k!
88         ssyr2k
89         dsyr2k
90         csyr2k
91         zsyr2k
92
93         unsafe-cher2k!
94         unsafe-zher2k!
95         cher2k!
96         zher2k!
97         cher2k
98         zher2k
99
100         unsafe-strmm!
101         unsafe-dtrmm!
102         unsafe-ctrmm!
103         unsafe-ztrmm!
104         strmm!
105         dtrmm!
106         ctrmm!
107         ztrmm!
108         strmm
109         dtrmm
110         ctrmm
111         ztrmm
112
113         unsafe-strsm!
114         unsafe-dtrsm!
115         unsafe-ctrsm!
116         unsafe-ztrsm!
117         strsm!
118         dtrsm!
119         ctrsm!
120         ztrsm!
121         strsm
122         dtrsm
123         ctrsm
124         ztrsm
125
126         unsafe-sgemv!
127         unsafe-dgemv!
128         unsafe-cgemv!
129         unsafe-zgemv!
130         sgemv!
131         dgemv!
132         cgemv!
133         zgemv!
134         sgemv
135         dgemv
136         cgemv
137         zgemv
138
139         unsafe-chemv!
140         unsafe-zhemv!
141         chemv!
142         zhemv!
143         chemv
144         zhemv
145         unsafe-chbmv!
146         unsafe-zhbmv!
147         chbmv!
148         zhbmv!
149         chbmv
150         zhbmv
151         unsafe-chpmv!
152         unsafe-zhpmv!
153         chpmv!
154         zhpmv!
155         chpmv
156         zhpmv
157         unsafe-ssymv!
158         unsafe-dsymv!
159         ssymv!
160         dsymv!
161         ssymv
162         dsymv
163         unsafe-ssbmv!
164         unsafe-dsbmv!
165         ssbmv!
166         dsbmv!
167         ssbmv
168         dsbmv
169         unsafe-sspmv!
170         unsafe-dspmv!
171         sspmv!
172         dspmv!
173         sspmv
174         dspmv
175         unsafe-strmv!
176         unsafe-dtrmv!
177         unsafe-ctrmv!
178         unsafe-ztrmv!
179         strmv!
180         dtrmv!
181         ctrmv!
182         ztrmv!
183         strmv
184         dtrmv
185         ctrmv
186         ztrmv
187         unsafe-stbmv!
188         unsafe-dtbmv!
189         unsafe-ctbmv!
190         unsafe-ztbmv!
191         stbmv!
192         dtbmv!
193         ctbmv!
194         ztbmv!
195         stbmv
196         dtbmv
197         ctbmv
198         ztbmv
199         unsafe-stpmv!
200         unsafe-dtpmv!
201         unsafe-ctpmv!
202         unsafe-ztpmv!
203         stpmv!
204         dtpmv!
205         ctpmv!
206         ztpmv!
207         stpmv
208         dtpmv
209         ctpmv
210         ztpmv
211         unsafe-strsv!
212         unsafe-dtrsv!
213         unsafe-ctrsv!
214         unsafe-ztrsv!
215         strsv!
216         dtrsv!
217         ctrsv!
218         ztrsv!
219         strsv
220         dtrsv
221         ctrsv
222         ztrsv
223         unsafe-stbsv!
224         unsafe-dtbsv!
225         unsafe-ctbsv!
226         unsafe-ztbsv!
227         stbsv!
228         dtbsv!
229         ctbsv!
230         ztbsv!
231         stbsv
232         dtbsv
233         ctbsv
234         ztbsv
235         unsafe-stpsv!
236         unsafe-dtpsv!
237         unsafe-ctpsv!
238         unsafe-ztpsv!
239         stpsv!
240         dtpsv!
241         ctpsv!
242         ztpsv!
243         stpsv
244         dtpsv
245         ctpsv
246         ztpsv
247         unsafe-sger!
248         unsafe-dger!
249         sger!
250         dger!
251         sger
252         dger
253         unsafe-siger!
254         unsafe-diger!
255         siger!
256         diger!
257         siger
258         diger
259         unsafe-cgeru!
260         unsafe-zgeru!
261         cgeru!
262         zgeru!
263         cgeru
264         zgeru
265         unsafe-cgerc!
266         unsafe-zgerc!
267         cgerc!
268         zgerc!
269         cgerc
270         zgerc
271         unsafe-cher!
272         unsafe-zher!
273         cher!
274         zher!
275         cher
276         zher
277         unsafe-chpr!
278         unsafe-zhpr!
279         chpr!
280         zhpr!
281         chpr
282         zhpr
283         unsafe-cher2!
284         unsafe-zher2!
285         cher2!
286         zher2!
287         cher2
288         zher2
289         unsafe-chpr2!
290         unsafe-zhpr2!
291         chpr2!
292         zhpr2!
293         chpr2
294         zhpr2
295         unsafe-ssyr!
296         unsafe-dsyr!
297         ssyr!
298         dsyr!
299         ssyr
300         dsyr
301         unsafe-sspr!
302         unsafe-dspr!
303         sspr!
304         dspr!
305         sspr
306         dspr
307         unsafe-ssyr2!
308         unsafe-dsyr2!
309         ssyr2!
310         dsyr2!
311         ssyr2
312         dsyr2
313         unsafe-sspr2!
314         unsafe-dspr2!
315         sspr2!
316         dspr2!
317         sspr2
318         dspr2
319
320         unsafe-srot!
321         unsafe-drot!
322         srot!
323         drot!
324         srot
325         drot
326         unsafe-srotm!
327         unsafe-drotm!
328         srotm!
329         drotm!
330         srotm
331         drotm
332         unsafe-sswap!
333         unsafe-dswap!
334         unsafe-cswap!
335         unsafe-zswap!
336         sswap!
337         dswap!
338         cswap!
339         zswap!
340         sswap
341         dswap
342         cswap
343         zswap
344         unsafe-sscal!
345         unsafe-dscal!
346         unsafe-cscal!
347         unsafe-zscal!
348         sscal!
349         dscal!
350         cscal!
351         zscal!
352         sscal
353         dscal
354         cscal
355         zscal
356         unsafe-saxpy!
357         unsafe-daxpy!
358         unsafe-caxpy!
359         unsafe-zaxpy!
360         saxpy!
361         daxpy!
362         caxpy!
363         zaxpy!
364         saxpy
365         daxpy
366         caxpy
367         zaxpy
368         unsafe-siaxpy!
369         unsafe-diaxpy!
370         unsafe-ciaxpy!
371         unsafe-ziaxpy!
372         siaxpy!
373         diaxpy!
374         ciaxpy!
375         ziaxpy!
376         siaxpy
377         diaxpy
378         ciaxpy
379         ziaxpy
380         sdot
381         ddot
382         cdotu
383         zdotu
384         cdotc
385         zdotc
386         snrm2
387         dnrm2
388         cnrm2
389         znrm2
390         sasum
391         dasum
392         casum
393         zasum
394         samax
395         damax
396         camax
397         zamax
398         )
399
400  (import scheme chicken data-structures foreign)
401
402  (require-extension srfi-4 bind)
403
404(define (blas:error x . rest)
405  (let ((port (open-output-string)))
406    (let loop ((objs (if (symbol? x) rest (cons x rest))))
407      (if (null? objs)
408          (begin
409            (newline port)
410            (error (if (symbol? x) x 'blas) 
411                   (get-output-string port)))
412          (begin (display (car objs) port)
413                 (display " " port)
414                 (loop (cdr objs)))))))
415
416
417(bind* #<<EOF
418
419typedef float CCOMPLEX;
420typedef double ZCOMPLEX;
421typedef int CBLAS_INDEX;
422
423
424/*
425 * Enumerated and derived types
426 */
427
428enum CBLAS_ORDER {CblasRowMajor=101, CblasColMajor=102};
429enum CBLAS_TRANSPOSE {CblasNoTrans=111, CblasTrans=112, CblasConjTrans=113};
430enum CBLAS_UPLO {CblasUpper=121, CblasLower=122};
431enum CBLAS_DIAG {CblasNonUnit=131, CblasUnit=132};
432enum CBLAS_SIDE {CblasLeft=141, CblasRight=142};
433
434/*
435 * ===========================================================================
436 * Prototypes for level 1 BLAS routines
437 * ===========================================================================
438 */
439
440/* 
441 * Routines with standard 4 prefixes (s, d, c, z)
442 */
443void cblas_sswap(const int N, float *X, const int incX, 
444                 float *Y, const int incY);
445void cblas_scopy(const int N, const float *X, const int incX, 
446                 float *Y, const int incY);
447void cblas_saxpy(const int N, const float alpha, const float *X,
448                 const int incX, float *Y, const int incY);
449
450void cblas_dswap(const int N, double *X, const int incX, 
451                 double *Y, const int incY);
452void cblas_dcopy(const int N, const double *X, const int incX, 
453                 double *Y, const int incY);
454void cblas_daxpy(const int N, const double alpha, const double *X,
455                 const int incX, double *Y, const int incY);
456
457void cblas_cswap(const int N, CCOMPLEX *X, const int incX, 
458                 CCOMPLEX *Y, const int incY);
459void cblas_ccopy(const int N, const CCOMPLEX *X, const int incX, 
460                 CCOMPLEX *Y, const int incY);
461void cblas_caxpy(const int N, const CCOMPLEX *alpha, const CCOMPLEX *X,
462                 const int incX, CCOMPLEX *Y, const int incY);
463
464void cblas_zswap(const int N, ZCOMPLEX *X, const int incX, 
465                 ZCOMPLEX *Y, const int incY);
466void cblas_zcopy(const int N, const ZCOMPLEX *X, const int incX, 
467                 ZCOMPLEX *Y, const int incY);
468void cblas_zaxpy(const int N, const ZCOMPLEX *alpha, const ZCOMPLEX *X,
469                 const int incX, ZCOMPLEX *Y, const int incY);
470
471
472/* 
473 * Routines with S and D prefix only
474 */
475void cblas_srotg(float *a, float *b, float *c, float *s);
476void cblas_srotmg(float *d1, float *d2, float *b1, const float b2, float *P);
477void cblas_srot(const int N, float *X, const int incX,
478                float *Y, const int incY, const float c, const float s);
479void cblas_srotm(const int N, float *X, const int incX,
480                float *Y, const int incY, const float *P);
481
482void cblas_drotg(double *a, double *b, double *c, double *s);
483void cblas_drotmg(double *d1, double *d2, double *b1, const double b2, double *P);
484void cblas_drot(const int N, double *X, const int incX,
485                double *Y, const int incY, const double c, const double  s);
486void cblas_drotm(const int N, double *X, const int incX,
487                double *Y, const int incY, const double *P);
488
489
490/* 
491 * Routines with S D C Z CS and ZD prefixes
492 */
493void cblas_sscal(const int N, const float alpha, float *X, const int incX);
494void cblas_dscal(const int N, const double alpha, double *X, const int incX);
495void cblas_cscal(const int N, const CCOMPLEX *alpha, CCOMPLEX *X, const int incX);
496void cblas_zscal(const int N, const ZCOMPLEX *alpha, ZCOMPLEX *X, const int incX);
497void cblas_csscal(const int N, const float alpha, CCOMPLEX *X, const int incX);
498void cblas_zdscal(const int N, const double alpha, ZCOMPLEX *X, const int incX);
499
500
501/* Offset variations of the copy, axpy routines */
502
503void sicopy(const int N, const float *X, const int incX, const
504                  int offsetX, float *Y, const int incY, const int offsetY)
505{
506  cblas_scopy (N, X+offsetX, incX, Y+offsetY, incY);
507}
508
509void dicopy(const int N, const double *X, const int incX, const
510                  int offsetX, double *Y, const int incY, const int offsetY)
511{
512  cblas_dcopy (N, X+offsetX, incX, Y+offsetY, incY);
513}
514
515
516void cicopy(const int N, const CCOMPLEX *X, const int incX, const
517                  int offsetX, CCOMPLEX *Y, const int incY, const int offsetY)
518{
519  cblas_ccopy (N, X+(2*offsetX), incX, Y+(2*offsetY), incY);
520}
521
522void zicopy(const int N, const ZCOMPLEX *X, const int incX, const
523                  int offsetX, ZCOMPLEX *Y, const int incY, const int offsetY)
524{
525  cblas_zcopy (N, X+(2*offsetX), incX, Y+(2*offsetY), incY);
526}
527
528
529void cblas_siaxpy(const int N, const float alpha, 
530                        const float *X, const int incX, const int offsetX, 
531                        float *Y, const int incY, const int offsetY)
532{
533
534 cblas_saxpy(N, alpha, X+offsetX, incX, Y+offsetY, incY);
535}
536
537
538void cblas_diaxpy(const int N, const double alpha, 
539                        const double *X, const int incX, const int offsetX, 
540                        double *Y, const int incY, const int offsetY)
541{
542
543 cblas_daxpy(N, alpha, X+offsetX, incX, Y+offsetY, incY);
544}
545
546
547void cblas_ciaxpy(const int N, const CCOMPLEX *alpha, 
548                        const CCOMPLEX *X, const int incX, const int offsetX, 
549                        CCOMPLEX *Y, const int incY, const int offsetY)
550{
551 cblas_caxpy(N, alpha, X+(2*offsetX), incX, Y+(2*offsetY), incY);
552}
553
554
555void cblas_ziaxpy(const int N, const ZCOMPLEX *alpha, 
556                        const ZCOMPLEX *X, const int incX, const int offsetX, 
557                        ZCOMPLEX *Y, const int incY, const int offsetY)
558{
559 cblas_zaxpy(N, alpha, X+(2*offsetX), incX, Y+(2*offsetY), incY);
560}
561EOF
562)
563
564
565(define RowMajor  CblasRowMajor)
566(define ColMajor  CblasColMajor)
567(define NoTrans   CblasNoTrans)
568(define Trans     CblasTrans)
569(define ConjTrans CblasConjTrans)
570(define Upper     CblasUpper)
571(define Lower     CblasLower)
572(define NonUnit   CblasNonUnit)
573(define Unit      CblasUnit)
574(define Left      CblasLeft)
575(define Right     CblasRight)
576
577
578(define (scopy x)
579  (let ((n (f32vector-length x)))
580    (let ((y  (make-f32vector n)))
581      (cblas_scopy n x 1 y 1)
582      y)))
583
584(define (dcopy x)
585  (let ((n (f64vector-length x)))
586    (let ((y  (make-f64vector n)))
587      (cblas_dcopy n x 1 y 1)
588      y)))
589
590(define (ccopy x)
591  (let ((n (fx/ (f32vector-length x) 2)))
592    (let ((y  (make-f32vector (fx* 2 n))))
593      (cblas_ccopy n x 1 y 1)
594      y)))
595
596(define (zcopy x)
597  (let ((n (fx/ (f64vector-length x) 2)))
598    (let ((y  (make-f64vector (fx* 2 n))))
599      (cblas_zcopy n x 1 y 1)
600      y)))
601
602
603(define-syntax icopy-wrapper
604  (lambda (x r c)
605    (let* ((copy            (cadr x))
606           (vector-length   (caddr x))
607           (make-vector     (cadddr x))
608           (name            copy)
609           (%define         (r 'define))
610           (%let            (r 'let))
611           (%cond           (r 'cond))
612           (%or             (r 'or))
613           (%if             (r 'if))
614           (%let-optionals  (r 'let-optionals)))
615         
616    `(,%define (,name n x . rest)
617       (,%let-optionals rest ((y #f) (offsetX 0) (offsetY 0) (incX 1) (incY 1))
618         (,%let ((xlen  (,vector-length x))
619               (ylen (,%if y (,vector-length y) (fx- n offsetX))))
620             (,%cond ((not (fx= n xlen))
621                    (blas:error ',name " n is not equal to the length of X (" xlen ")"))
622                   ((fx< offsetX 0) 
623                    (blas:error ',name "offset of vector X (" offsetX ") is negative"))
624                   ((fx>= offsetX xlen)
625                    (blas:error ',name "offset of vector X (" offsetX ") is greater than or equal to its length: " xlen))
626                   ((fx< offsetX 0) 
627                    (blas:error ',name "offset of vector X (" offsetX ") is negative"))
628                   ((fx>= offsetY ylen)
629                    (blas:error ',name "offset of vector Y (" offsetY ") is greater than or equal to its length: " ylen))
630                   ((fx> (- ylen offsetY) (- xlen offsetX))
631                    (blas:error ',name "range of vector Y (" (- ylen offsetY) 
632                                ") is greater than range of vector X: " ( - xlen offsetX))))
633             (,%let ((y (,%or y (,make-vector ylen))))
634               (,copy n x incX offsetX y incY offsetY)
635               y))))))
636  )
637
638(icopy-wrapper sicopy f32vector-length make-f32vector)
639(icopy-wrapper dicopy f64vector-length make-f64vector)
640(icopy-wrapper cicopy 
641               (lambda (x) (fx/ (f32vector-length x) 2))
642               (lambda (n) (make-f32vector (fx* 2 n))))
643(icopy-wrapper zicopy 
644               (lambda (x) (fx/ (f64vector-length x) 2))
645               (lambda (n) (make-f64vector (fx* 2 n))))
646
647
648
649
650
651(bind* #<<EOF
652
653/*
654 * ===========================================================================
655 * Prototypes for level 1 BLAS functions (complex are recast as routines)
656 * ===========================================================================
657 */
658float  cblas_sdsdot(const int N, const float alpha, const float *X,
659                    const int incX, const float *Y, const int incY);
660double cblas_dsdot(const int N, const float *X, const int incX, const float *Y,
661                   const int incY);
662float  cblas_sdot(const int N, const float  *X, const int incX,
663                  const float  *Y, const int incY);
664double cblas_ddot(const int N, const double *X, const int incX,
665                  const double *Y, const int incY);
666
667/*
668 * Functions having prefixes Z and C only
669 */
670void   cblas_cdotu_sub(const int N, const CCOMPLEX *X, const int incX,
671                       const CCOMPLEX *Y, const int incY, CCOMPLEX *dotu);
672void   cblas_cdotc_sub(const int N, const CCOMPLEX *X, const int incX,
673                       const CCOMPLEX *Y, const int incY, CCOMPLEX *dotc);
674
675void   cblas_zdotu_sub(const int N, const ZCOMPLEX *X, const int incX,
676                       const ZCOMPLEX *Y, const int incY, ZCOMPLEX *dotu);
677void   cblas_zdotc_sub(const int N, const ZCOMPLEX *X, const int incX,
678                       const ZCOMPLEX *Y, const int incY, ZCOMPLEX *dotc);
679
680
681/*
682 * Functions having prefixes S D SC DZ
683 */
684float  cblas_snrm2(const int N, const float *X, const int incX);
685float  cblas_sasum(const int N, const float *X, const int incX);
686
687double cblas_dnrm2(const int N, const double *X, const int incX);
688double cblas_dasum(const int N, const double *X, const int incX);
689
690float  cblas_scnrm2(const int N, const CCOMPLEX *X, const int incX);
691float  cblas_scasum(const int N, const CCOMPLEX *X, const int incX);
692
693double cblas_dznrm2(const int N, const ZCOMPLEX *X, const int incX);
694double cblas_dzasum(const int N, const ZCOMPLEX *X, const int incX);
695
696
697
698/*
699 * Functions having standard 4 prefixes (S D C Z)
700 */
701CBLAS_INDEX cblas_isamax(const int N, const float  *X, const int incX);
702CBLAS_INDEX cblas_idamax(const int N, const double *X, const int incX);
703CBLAS_INDEX cblas_icamax(const int N, const void   *X, const int incX);
704CBLAS_INDEX cblas_izamax(const int N, const void   *X, const int incX);
705
706/*
707 * ===========================================================================
708 * Prototypes for level 2 BLAS
709 * ===========================================================================
710 */
711
712/* 
713 * Routines with standard 4 prefixes (S, D, C, Z)
714 */
715void cblas_sgemv(const enum CBLAS_ORDER order,
716                 const enum CBLAS_TRANSPOSE TransA, const int M, const int N,
717                 const float alpha, const float *A, const int lda,
718                 const float *X, const int incX, const float beta,
719                 float *Y, const int incY);
720void cblas_sgbmv(const enum CBLAS_ORDER order,
721                 const enum CBLAS_TRANSPOSE TransA, const int M, const int N,
722                 const int KL, const int KU, const float alpha,
723                 const float *A, const int lda, const float *X,
724                 const int incX, const float beta, float *Y, const int incY);
725void cblas_strmv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
726                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
727                 const int N, const float *A, const int lda, 
728                 float *X, const int incX);
729void cblas_stbmv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
730                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
731                 const int N, const int K, const float *A, const int lda, 
732                 float *X, const int incX);
733void cblas_stpmv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
734                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
735                 const int N, const float *Ap, float *X, const int incX);
736void cblas_strsv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
737                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
738                 const int N, const float *A, const int lda, float *X,
739                 const int incX);
740void cblas_stbsv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
741                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
742                 const int N, const int K, const float *A, const int lda,
743                 float *X, const int incX);
744void cblas_stpsv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
745                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
746                 const int N, const float *Ap, float *X, const int incX);
747
748void cblas_dgemv(const enum CBLAS_ORDER order,
749                 const enum CBLAS_TRANSPOSE TransA, const int M, const int N,
750                 const double alpha, const double *A, const int lda,
751                 const double *X, const int incX, const double beta,
752                 double *Y, const int incY);
753void cblas_dgbmv(const enum CBLAS_ORDER order,
754                 const enum CBLAS_TRANSPOSE TransA, const int M, const int N,
755                 const int KL, const int KU, const double alpha,
756                 const double *A, const int lda, const double *X,
757                 const int incX, const double beta, double *Y, const int incY);
758void cblas_dtrmv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
759                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
760                 const int N, const double *A, const int lda, 
761                 double *X, const int incX);
762void cblas_dtbmv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
763                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
764                 const int N, const int K, const double *A, const int lda, 
765                 double *X, const int incX);
766void cblas_dtpmv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
767                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
768                 const int N, const double *Ap, double *X, const int incX);
769void cblas_dtrsv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
770                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
771                 const int N, const double *A, const int lda, double *X,
772                 const int incX);
773void cblas_dtbsv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
774                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
775                 const int N, const int K, const double *A, const int lda,
776                 double *X, const int incX);
777void cblas_dtpsv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
778                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
779                 const int N, const double *Ap, double *X, const int incX);
780
781void cblas_cgemv(const enum CBLAS_ORDER order,
782                 const enum CBLAS_TRANSPOSE TransA, const int M, const int N,
783                 const CCOMPLEX *alpha, const CCOMPLEX *A, const int lda,
784                 const CCOMPLEX *X, const int incX, const CCOMPLEX *beta,
785                 CCOMPLEX *Y, const int incY);
786void cblas_cgbmv(const enum CBLAS_ORDER order,
787                 const enum CBLAS_TRANSPOSE TransA, const int M, const int N,
788                 const int KL, const int KU, const CCOMPLEX *alpha,
789                 const CCOMPLEX *A, const int lda, const CCOMPLEX *X,
790                 const int incX, const CCOMPLEX *beta, CCOMPLEX *Y, const int incY);
791void cblas_ctrmv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
792                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
793                 const int N, const CCOMPLEX *A, const int lda, 
794                 CCOMPLEX *X, const int incX);
795void cblas_ctbmv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
796                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
797                 const int N, const int K, const CCOMPLEX *A, const int lda, 
798                 CCOMPLEX *X, const int incX);
799void cblas_ctpmv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
800                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
801                 const int N, const CCOMPLEX *Ap, CCOMPLEX *X, const int incX);
802void cblas_ctrsv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
803                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
804                 const int N, const CCOMPLEX *A, const int lda, CCOMPLEX *X,
805                 const int incX);
806void cblas_ctbsv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
807                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
808                 const int N, const int K, const CCOMPLEX *A, const int lda,
809                 CCOMPLEX *X, const int incX);
810void cblas_ctpsv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
811                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
812                 const int N, const CCOMPLEX *Ap, CCOMPLEX *X, const int incX);
813
814void cblas_zgemv(const enum CBLAS_ORDER order,
815                 const enum CBLAS_TRANSPOSE TransA, const int M, const int N,
816                 const ZCOMPLEX *alpha, const ZCOMPLEX *A, const int lda,
817                 const ZCOMPLEX *X, const int incX, const ZCOMPLEX *beta,
818                 ZCOMPLEX *Y, const int incY);
819void cblas_zgbmv(const enum CBLAS_ORDER order,
820                 const enum CBLAS_TRANSPOSE TransA, const int M, const int N,
821                 const int KL, const int KU, const ZCOMPLEX *alpha,
822                 const ZCOMPLEX *A, const int lda, const ZCOMPLEX *X,
823                 const int incX, const ZCOMPLEX *beta, ZCOMPLEX *Y, const int incY);
824void cblas_ztrmv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
825                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
826                 const int N, const ZCOMPLEX *A, const int lda, 
827                 ZCOMPLEX *X, const int incX);
828void cblas_ztbmv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
829                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
830                 const int N, const int K, const ZCOMPLEX *A, const int lda, 
831                 ZCOMPLEX *X, const int incX);
832void cblas_ztpmv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
833                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
834                 const int N, const ZCOMPLEX *Ap, ZCOMPLEX *X, const int incX);
835void cblas_ztrsv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
836                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
837                 const int N, const ZCOMPLEX *A, const int lda, ZCOMPLEX *X,
838                 const int incX);
839void cblas_ztbsv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
840                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
841                 const int N, const int K, const ZCOMPLEX *A, const int lda,
842                 ZCOMPLEX *X, const int incX);
843void cblas_ztpsv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
844                 const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_DIAG Diag,
845                 const int N, const ZCOMPLEX *Ap, ZCOMPLEX *X, const int incX);
846
847
848/* 
849 * Routines with S and D prefixes only
850 */
851void cblas_ssymv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
852                 const int N, const float alpha, const float *A,
853                 const int lda, const float *X, const int incX,
854                 const float beta, float *Y, const int incY);
855void cblas_ssbmv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
856                 const int N, const int K, const float alpha, const float *A,
857                 const int lda, const float *X, const int incX,
858                 const float beta, float *Y, const int incY);
859void cblas_sspmv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
860                 const int N, const float alpha, const float *Ap,
861                 const float *X, const int incX,
862                 const float beta, float *Y, const int incY);
863void cblas_sger(const enum CBLAS_ORDER order, const int M, const int N,
864                const float alpha, const float *X, const int incX,
865                const float *Y, const int incY, float *A, const int lda);
866void cblas_ssyr(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
867                const int N, const float alpha, const float *X,
868                const int incX, float *A, const int lda);
869void cblas_sspr(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
870                const int N, const float alpha, const float *X,
871                const int incX, float *Ap);
872void cblas_ssyr2(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
873                const int N, const float alpha, const float *X,
874                const int incX, const float *Y, const int incY, float *A,
875                const int lda);
876void cblas_sspr2(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
877                const int N, const float alpha, const float *X,
878                const int incX, const float *Y, const int incY, float *A);
879
880void cblas_dsymv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
881                 const int N, const double alpha, const double *A,
882                 const int lda, const double *X, const int incX,
883                 const double beta, double *Y, const int incY);
884void cblas_dsbmv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
885                 const int N, const int K, const double alpha, const double *A,
886                 const int lda, const double *X, const int incX,
887                 const double beta, double *Y, const int incY);
888void cblas_dspmv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
889                 const int N, const double alpha, const double *Ap,
890                 const double *X, const int incX,
891                 const double beta, double *Y, const int incY);
892void cblas_dger(const enum CBLAS_ORDER order, const int M, const int N,
893                const double alpha, const double *X, const int incX,
894                const double *Y, const int incY, double *A, const int lda);
895void cblas_dsyr(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
896                const int N, const double alpha, const double *X,
897                const int incX, double *A, const int lda);
898void cblas_dspr(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
899                const int N, const double alpha, const double *X,
900                const int incX, double *Ap);
901void cblas_dsyr2(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
902                const int N, const double alpha, const double *X,
903                const int incX, const double *Y, const int incY, double *A,
904                const int lda);
905void cblas_dspr2(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
906                const int N, const double alpha, const double *X,
907                const int incX, const double *Y, const int incY, double *A);
908
909
910/* 
911 * Routines with C and Z prefixes only
912 */
913void cblas_chemv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
914                 const int N, const CCOMPLEX *alpha, const CCOMPLEX *A,
915                 const int lda, const CCOMPLEX *X, const int incX,
916                 const CCOMPLEX *beta, CCOMPLEX *Y, const int incY);
917void cblas_chbmv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
918                 const int N, const int K, const CCOMPLEX *alpha, const CCOMPLEX *A,
919                 const int lda, const CCOMPLEX *X, const int incX,
920                 const CCOMPLEX *beta, CCOMPLEX *Y, const int incY);
921void cblas_chpmv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
922                 const int N, const CCOMPLEX *alpha, const CCOMPLEX *Ap,
923                 const CCOMPLEX *X, const int incX,
924                 const CCOMPLEX *beta, CCOMPLEX *Y, const int incY);
925void cblas_cgeru(const enum CBLAS_ORDER order, const int M, const int N,
926                 const CCOMPLEX *alpha, const CCOMPLEX *X, const int incX,
927                 const CCOMPLEX *Y, const int incY, CCOMPLEX *A, const int lda);
928void cblas_cgerc(const enum CBLAS_ORDER order, const int M, const int N,
929                 const CCOMPLEX *alpha, const CCOMPLEX *X, const int incX,
930                 const CCOMPLEX *Y, const int incY, CCOMPLEX *A, const int lda);
931void cblas_cher(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
932                const int N, const float alpha, const CCOMPLEX *X, const int incX,
933                CCOMPLEX *A, const int lda);
934void cblas_chpr(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
935                const int N, const float alpha, const CCOMPLEX *X,
936                const int incX, CCOMPLEX *A);
937void cblas_cher2(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo, const int N,
938                const CCOMPLEX *alpha, const CCOMPLEX *X, const int incX,
939                const CCOMPLEX *Y, const int incY, CCOMPLEX *A, const int lda);
940void cblas_chpr2(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo, const int N,
941                const CCOMPLEX *alpha, const CCOMPLEX *X, const int incX,
942                const CCOMPLEX *Y, const int incY, CCOMPLEX *Ap);
943
944void cblas_zhemv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
945                 const int N, const ZCOMPLEX *alpha, const ZCOMPLEX *A,
946                 const int lda, const ZCOMPLEX *X, const int incX,
947                 const ZCOMPLEX *beta, ZCOMPLEX *Y, const int incY);
948void cblas_zhbmv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
949                 const int N, const int K, const ZCOMPLEX *alpha, const ZCOMPLEX *A,
950                 const int lda, const ZCOMPLEX *X, const int incX,
951                 const ZCOMPLEX *beta, ZCOMPLEX *Y, const int incY);
952void cblas_zhpmv(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
953                 const int N, const ZCOMPLEX *alpha, const ZCOMPLEX *Ap,
954                 const ZCOMPLEX *X, const int incX,
955                 const ZCOMPLEX *beta, ZCOMPLEX *Y, const int incY);
956void cblas_zgeru(const enum CBLAS_ORDER order, const int M, const int N,
957                 const ZCOMPLEX *alpha, const ZCOMPLEX *X, const int incX,
958                 const ZCOMPLEX *Y, const int incY, ZCOMPLEX *A, const int lda);
959void cblas_zgerc(const enum CBLAS_ORDER order, const int M, const int N,
960                 const ZCOMPLEX *alpha, const ZCOMPLEX *X, const int incX,
961                 const ZCOMPLEX *Y, const int incY, ZCOMPLEX *A, const int lda);
962void cblas_zher(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
963                const int N, const double alpha, const ZCOMPLEX *X, const int incX,
964                ZCOMPLEX *A, const int lda);
965void cblas_zhpr(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
966                const int N, const double alpha, const ZCOMPLEX *X,
967                const int incX, ZCOMPLEX *A);
968void cblas_zher2(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo, const int N,
969                const ZCOMPLEX *alpha, const ZCOMPLEX *X, const int incX,
970                const ZCOMPLEX *Y, const int incY, ZCOMPLEX *A, const int lda);
971void cblas_zhpr2(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo, const int N,
972                const ZCOMPLEX *alpha, const ZCOMPLEX *X, const int incX,
973                const ZCOMPLEX *Y, const int incY, ZCOMPLEX *Ap);
974
975/*
976 * ===========================================================================
977 * Prototypes for level 3 BLAS
978 * ===========================================================================
979 */
980
981/* 
982 * Routines with standard 4 prefixes (S, D, C, Z)
983 */
984void cblas_sgemm(const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA,
985                 const enum CBLAS_TRANSPOSE TransB, const int M, const int N,
986                 const int K, const float alpha, const float *A,
987                 const int lda, const float *B, const int ldb,
988                 const float beta, float *C, const int ldc);
989void cblas_ssymm(const enum CBLAS_ORDER Order, const enum CBLAS_SIDE Side,
990                 const enum CBLAS_UPLO Uplo, const int M, const int N,
991                 const float alpha, const float *A, const int lda,
992                 const float *B, const int ldb, const float beta,
993                 float *C, const int ldc);
994void cblas_ssyrk(const enum CBLAS_ORDER Order, const enum CBLAS_UPLO Uplo,
995                 const enum CBLAS_TRANSPOSE Trans, const int N, const int K,
996                 const float alpha, const float *A, const int lda,
997                 const float beta, float *C, const int ldc);
998void cblas_ssyr2k(const enum CBLAS_ORDER Order, const enum CBLAS_UPLO Uplo,
999                  const enum CBLAS_TRANSPOSE Trans, const int N, const int K,
1000                  const float alpha, const float *A, const int lda,
1001                  const float *B, const int ldb, const float beta,
1002                  float *C, const int ldc);
1003void cblas_strmm(const enum CBLAS_ORDER Order, const enum CBLAS_SIDE Side,
1004                 const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA,
1005                 const enum CBLAS_DIAG Diag, const int M, const int N,
1006                 const float alpha, const float *A, const int lda,
1007                 float *B, const int ldb);
1008void cblas_strsm(const enum CBLAS_ORDER Order, const enum CBLAS_SIDE Side,
1009                 const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA,
1010                 const enum CBLAS_DIAG Diag, const int M, const int N,
1011                 const float alpha, const float *A, const int lda,
1012                 float *B, const int ldb);
1013
1014void cblas_dgemm(const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA,
1015                 const enum CBLAS_TRANSPOSE TransB, const int M, const int N,
1016                 const int K, const double alpha, const double *A,
1017                 const int lda, const double *B, const int ldb,
1018                 const double beta, double *C, const int ldc);
1019void cblas_dsymm(const enum CBLAS_ORDER Order, const enum CBLAS_SIDE Side,
1020                 const enum CBLAS_UPLO Uplo, const int M, const int N,
1021                 const double alpha, const double *A, const int lda,
1022                 const double *B, const int ldb, const double beta,
1023                 double *C, const int ldc);
1024void cblas_dsyrk(const enum CBLAS_ORDER Order, const enum CBLAS_UPLO Uplo,
1025                 const enum CBLAS_TRANSPOSE Trans, const int N, const int K,
1026                 const double alpha, const double *A, const int lda,
1027                 const double beta, double *C, const int ldc);
1028void cblas_dsyr2k(const enum CBLAS_ORDER Order, const enum CBLAS_UPLO Uplo,
1029                  const enum CBLAS_TRANSPOSE Trans, const int N, const int K,
1030                  const double alpha, const double *A, const int lda,
1031                  const double *B, const int ldb, const double beta,
1032                  double *C, const int ldc);
1033void cblas_dtrmm(const enum CBLAS_ORDER Order, const enum CBLAS_SIDE Side,
1034                 const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA,
1035                 const enum CBLAS_DIAG Diag, const int M, const int N,
1036                 const double alpha, const double *A, const int lda,
1037                 double *B, const int ldb);
1038void cblas_dtrsm(const enum CBLAS_ORDER Order, const enum CBLAS_SIDE Side,
1039                 const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA,
1040                 const enum CBLAS_DIAG Diag, const int M, const int N,
1041                 const double alpha, const double *A, const int lda,
1042                 double *B, const int ldb);
1043
1044void cblas_cgemm(const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA,
1045                 const enum CBLAS_TRANSPOSE TransB, const int M, const int N,
1046                 const int K, const CCOMPLEX *alpha, const CCOMPLEX *A,
1047                 const int lda, const CCOMPLEX *B, const int ldb,
1048                 const CCOMPLEX *beta, CCOMPLEX *C, const int ldc);
1049void cblas_csymm(const enum CBLAS_ORDER Order, const enum CBLAS_SIDE Side,
1050                 const enum CBLAS_UPLO Uplo, const int M, const int N,
1051                 const CCOMPLEX *alpha, const CCOMPLEX *A, const int lda,
1052                 const CCOMPLEX *B, const int ldb, const CCOMPLEX *beta,
1053                 CCOMPLEX *C, const int ldc);
1054void cblas_csyrk(const enum CBLAS_ORDER Order, const enum CBLAS_UPLO Uplo,
1055                 const enum CBLAS_TRANSPOSE Trans, const int N, const int K,
1056                 const CCOMPLEX *alpha, const CCOMPLEX *A, const int lda,
1057                 const CCOMPLEX *beta, CCOMPLEX *C, const int ldc);
1058void cblas_csyr2k(const enum CBLAS_ORDER Order, const enum CBLAS_UPLO Uplo,
1059                  const enum CBLAS_TRANSPOSE Trans, const int N, const int K,
1060                  const CCOMPLEX *alpha, const CCOMPLEX *A, const int lda,
1061                  const CCOMPLEX *B, const int ldb, const CCOMPLEX *beta,
1062                  CCOMPLEX *C, const int ldc);
1063void cblas_ctrmm(const enum CBLAS_ORDER Order, const enum CBLAS_SIDE Side,
1064                 const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA,
1065                 const enum CBLAS_DIAG Diag, const int M, const int N,
1066                 const CCOMPLEX *alpha, const CCOMPLEX *A, const int lda,
1067                 CCOMPLEX *B, const int ldb);
1068void cblas_ctrsm(const enum CBLAS_ORDER Order, const enum CBLAS_SIDE Side,
1069                 const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA,
1070                 const enum CBLAS_DIAG Diag, const int M, const int N,
1071                 const CCOMPLEX *alpha, const CCOMPLEX *A, const int lda,
1072                 CCOMPLEX *B, const int ldb);
1073
1074void cblas_zgemm(const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA,
1075                 const enum CBLAS_TRANSPOSE TransB, const int M, const int N,
1076                 const int K, const ZCOMPLEX *alpha, const ZCOMPLEX *A,
1077                 const int lda, const ZCOMPLEX *B, const int ldb,
1078                 const ZCOMPLEX *beta, ZCOMPLEX *C, const int ldc);
1079void cblas_zsymm(const enum CBLAS_ORDER Order, const enum CBLAS_SIDE Side,
1080                 const enum CBLAS_UPLO Uplo, const int M, const int N,
1081                 const ZCOMPLEX *alpha, const ZCOMPLEX *A, const int lda,
1082                 const ZCOMPLEX *B, const int ldb, const ZCOMPLEX *beta,
1083                 ZCOMPLEX *C, const int ldc);
1084void cblas_zsyrk(const enum CBLAS_ORDER Order, const enum CBLAS_UPLO Uplo,
1085                 const enum CBLAS_TRANSPOSE Trans, const int N, const int K,
1086                 const ZCOMPLEX *alpha, const ZCOMPLEX *A, const int lda,
1087                 const ZCOMPLEX *beta, ZCOMPLEX *C, const int ldc);
1088void cblas_zsyr2k(const enum CBLAS_ORDER Order, const enum CBLAS_UPLO Uplo,
1089                  const enum CBLAS_TRANSPOSE Trans, const int N, const int K,
1090                  const ZCOMPLEX *alpha, const ZCOMPLEX *A, const int lda,
1091                  const ZCOMPLEX *B, const int ldb, const ZCOMPLEX *beta,
1092                  ZCOMPLEX *C, const int ldc);
1093void cblas_ztrmm(const enum CBLAS_ORDER Order, const enum CBLAS_SIDE Side,
1094                 const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA,
1095                 const enum CBLAS_DIAG Diag, const int M, const int N,
1096                 const ZCOMPLEX *alpha, const ZCOMPLEX *A, const int lda,
1097                 ZCOMPLEX *B, const int ldb);
1098void cblas_ztrsm(const enum CBLAS_ORDER Order, const enum CBLAS_SIDE Side,
1099                 const enum CBLAS_UPLO Uplo, const enum CBLAS_TRANSPOSE TransA,
1100                 const enum CBLAS_DIAG Diag, const int M, const int N,
1101                 const ZCOMPLEX *alpha, const ZCOMPLEX *A, const int lda,
1102                 ZCOMPLEX *B, const int ldb);
1103
1104
1105/* 
1106 * Routines with prefixes C and Z only
1107 */
1108void cblas_chemm(const enum CBLAS_ORDER Order, const enum CBLAS_SIDE Side,
1109                 const enum CBLAS_UPLO Uplo, const int M, const int N,
1110                 const CCOMPLEX *alpha, const CCOMPLEX *A, const int lda,
1111                 const CCOMPLEX *B, const int ldb, const CCOMPLEX *beta,
1112                 CCOMPLEX *C, const int ldc);
1113void cblas_cherk(const enum CBLAS_ORDER Order, const enum CBLAS_UPLO Uplo,
1114                 const enum CBLAS_TRANSPOSE Trans, const int N, const int K,
1115                 const float alpha, const CCOMPLEX *A, const int lda,
1116                 const float beta, CCOMPLEX *C, const int ldc);
1117void cblas_cher2k(const enum CBLAS_ORDER Order, const enum CBLAS_UPLO Uplo,
1118                  const enum CBLAS_TRANSPOSE Trans, const int N, const int K,
1119                  const CCOMPLEX *alpha, const CCOMPLEX *A, const int lda,
1120                  const CCOMPLEX *B, const int ldb, const float beta,
1121                  CCOMPLEX *C, const int ldc);
1122
1123void cblas_zhemm(const enum CBLAS_ORDER Order, const enum CBLAS_SIDE Side,
1124                 const enum CBLAS_UPLO Uplo, const int M, const int N,
1125                 const ZCOMPLEX *alpha, const ZCOMPLEX *A, const int lda,
1126                 const ZCOMPLEX *B, const int ldb, const ZCOMPLEX *beta,
1127                 ZCOMPLEX *C, const int ldc);
1128void cblas_zherk(const enum CBLAS_ORDER Order, const enum CBLAS_UPLO Uplo,
1129                 const enum CBLAS_TRANSPOSE Trans, const int N, const int K,
1130                 const double alpha, const ZCOMPLEX *A, const int lda,
1131                 const double beta, ZCOMPLEX *C, const int ldc);
1132void cblas_zher2k(const enum CBLAS_ORDER Order, const enum CBLAS_UPLO Uplo,
1133                  const enum CBLAS_TRANSPOSE Trans, const int N, const int K,
1134                  const ZCOMPLEX *alpha, const ZCOMPLEX *A, const int lda,
1135                  const ZCOMPLEX *B, const int ldb, const double beta,
1136                  ZCOMPLEX *C, const int ldc);
1137
1138
1139/* Offset variants of ger routines */
1140
1141void cblas_siger(const enum CBLAS_ORDER order, const int M, const int N,
1142                       const float alpha, 
1143                       const float *X, const int incX, const int offsetX,
1144                       const float *Y, const int incY, const int offsetY,
1145                       float *A, const int lda)
1146{
1147
1148 cblas_sger(order, M, N, alpha, X+offsetX, incX, Y+offsetY, incY, A, lda);
1149}
1150
1151void cblas_diger(const enum CBLAS_ORDER order, const int M, const int N,
1152                       const double alpha, 
1153                       const double *X, const int incX, const int offsetX,
1154                       const double *Y, const int incY, const int offsetY,
1155                       double *A, const int lda)
1156{
1157
1158 cblas_dger(order, M, N, alpha, X+offsetX, incX, Y+offsetY, incY, A, lda);
1159}
1160
1161
1162/* "Standardize" some procedure names */
1163
1164float  cblas_cnrm2(const int N, const CCOMPLEX *X, const int incX)
1165{
1166 return cblas_scnrm2(N,X,incX);
1167}
1168
1169double cblas_znrm2(const int N, const ZCOMPLEX *X, const int incX)
1170{
1171 return cblas_dznrm2(N,X,incX);
1172}
1173
1174
1175float  cblas_casum(const int N, const CCOMPLEX *X, const int incX)
1176{
1177  return cblas_scasum(N,X,incX);
1178}
1179
1180double  cblas_zasum(const int N, const ZCOMPLEX *X, const int incX)
1181{
1182  return cblas_dzasum(N,X,incX);
1183}
1184
1185CBLAS_INDEX cblas_samax(const int N, const float  *X, const int incX)
1186{
1187 return cblas_isamax(N,X,incX);
1188}
1189
1190CBLAS_INDEX cblas_damax(const int N, const double *X, const int incX)
1191{
1192 return cblas_idamax(N,X,incX);
1193}
1194
1195CBLAS_INDEX cblas_camax(const int N, const void   *X, const int incX)
1196{
1197 return cblas_icamax(N,X,incX);
1198}
1199
1200CBLAS_INDEX cblas_zamax(const int N, const void   *X, const int incX)
1201{
1202 return cblas_izamax(N,X,incX);
1203}
1204
1205void   cblas_cdotu(const int N, const CCOMPLEX *X, const int incX,
1206                         const CCOMPLEX *Y, const int incY, CCOMPLEX *dotu)
1207{
1208 cblas_cdotu_sub(N,X,incX,Y,incY,dotu);
1209}
1210
1211void   cblas_cdotc(const int N, const CCOMPLEX *X, const int incX,
1212                         const CCOMPLEX *Y, const int incY, CCOMPLEX *dotc)
1213{
1214 cblas_cdotc_sub(N,X,incX,Y,incY,dotc);
1215}
1216
1217void   cblas_zdotu(const int N, const ZCOMPLEX *X, const int incX,
1218                       const ZCOMPLEX *Y, const int incY, ZCOMPLEX *dotu)
1219{
1220 cblas_zdotu_sub(N,X,incX,Y,incY,dotu);
1221}
1222
1223void   cblas_zdotc(const int N, const ZCOMPLEX *X, const int incX,
1224                       const ZCOMPLEX *Y, const int incY, ZCOMPLEX *dotc)
1225{
1226 cblas_zdotc_sub(N,X,incX,Y,incY,dotc);
1227}
1228EOF
1229)
1230
1231
1232
1233(define-syntax blas-level3-wrap
1234      (lambda (x r c)
1235        (let* ((fn      (cadr x))
1236               (ret     (caddr x))
1237               (err     (cadddr x))
1238               (vsize   (car (cddddr x)))
1239               (copy    (cadr (cddddr x)))
1240               (cfname  (string->symbol (conc "cblas_" (symbol->string (car fn)))))
1241               (fname   (string->symbol (conc (if vsize "" "unsafe-") 
1242                                              (symbol->string (car fn))
1243                                              (if copy "" "!"))))
1244               (%define         (r 'define))
1245               (%begin          (r 'begin))
1246               (%let            (r 'let))
1247               (%cond           (r 'cond))
1248               (%or             (r 'or))
1249               (%if             (r 'if))
1250               (%let-optionals  (r 'let-optionals))
1251               
1252               (ka              (r 'ka))
1253               (kb              (r 'kb))
1254               (kc              (r 'kc))
1255               (asize           (r 'asize))
1256               (bsize           (r 'bsize))
1257               (csize           (r 'csize))
1258               
1259               (args   (reverse (cdr fn)))
1260
1261               (fsig  (let loop ((args args) (sig 'rest))
1262                        (if (null? args) (cons fname sig)
1263                            (let ((x (car args)))
1264                              (let ((sig (case x 
1265                                          ((lda)     sig)
1266                                          ((ldb)     sig)
1267                                          ((ldc)     sig)
1268                                          (else      (cons x sig)))))
1269                                (loop (cdr args) sig))))))
1270
1271               (opts  (append
1272                       (if (memq 'lda fn) 
1273                           `((lda 
1274                              ,(cond ((memq 'side fn)
1275                                      `(,%if (= side Left) m n))
1276                                     ((memq 'transA fn)
1277                                      `(,%if (= transA NoTrans) k ,(if (memq 'm fn) 'm 'n)))
1278                                     ((memq 'trans fn)
1279                                      `(,%if (= trans NoTrans) k n))
1280                                     (else  
1281                                      (cond ((memq 'm fn) 'm) 
1282                                            (else 'n))))))
1283                           `())
1284                       (if (memq 'ldb fn)   
1285                           `((ldb ,(cond ((memq 'transB fn) 
1286                                          `(,%if (= transB NoTrans) n k))
1287                                         ((memq 'trans fn) 
1288                                          `(,%if (= trans NoTrans) k n))
1289                                         (else 'n))))
1290                           `())
1291                       (if (memq 'ldc fn) 
1292                           `((ldc  n)) `()))))
1293
1294          `(,%define ,fsig 
1295                     (,%let-optionals rest ,opts
1296                      ,(if vsize
1297                           `(,%begin
1298                             (,%let ((,asize (,vsize a))
1299                                     (,ka    ,(cond ((memq 'side fn) 
1300                                                     `(,%if (= side Left) m n))
1301                                                    ((memq 'transA fn) 
1302                                                     `(,%if (= transA NoTrans) 
1303                                                            ,(if (memq 'm fn) 'm 'n) k))
1304                                                    ((memq 'trans fn) 
1305                                                     `(,%if (= trans NoTrans) 
1306                                                            ,(if (memq 'm fn) 'm 'n) k))
1307                                                    (else (if (memq 'm fn) 'm 'n)))))
1308                                    (,%if (< ,asize (fx* lda ,ka)) 
1309                                          (blas:error ',fname (conc "matrix A is allocated " ,asize " elements "
1310                                                                    "but given dimensions are " ,ka " by " lda))))
1311                             ,(if (memq 'b fn)
1312                                  `(,%let ((,bsize (,vsize b))
1313                                           (,kb    ,(cond ((memq 'transB fn) 
1314                                                           `(,%if (= transB NoTrans) k n))
1315                                                          ((memq 'trans fn) 
1316                                                           `(,%if (= trans NoTrans) n k))
1317                                                          (else 'm))))
1318                                          (,%if (< ,bsize (fx* ldb ,kb)) 
1319                                                (blas:error ',fname (conc "matrix B is allocated " ,bsize " elements "
1320                                                                          "but given dimensions are " ,kb " by " ldb))))
1321                                  `(begin))
1322                             ,(if (memq 'c fn)
1323                                  `(let ((,csize (,vsize c))
1324                                         (,kc    ,(if (memq 'm fn) 'm 'n)))
1325                                     (if (< ,csize (fx* ldc ,kc))
1326                                         (blas:error ',fname (conc "matrix C is allocated " ,csize " elements "
1327                                                                   "but given dimensions are " ,kc " by " ldc))))
1328                                  `(begin)))
1329                           `(begin))
1330                      (,%let ,(let loop ((fn fn) (bnds '()))
1331                                (if (null? fn) bnds
1332                                    (let ((x (car fn)))
1333                                      (let ((bnds (case x 
1334                                                    (else    (if (and copy (memq x ret))
1335                                                                 (cons `(,x (,copy ,x)) bnds)
1336                                                                 bnds)))))
1337                                        (loop (cdr fn) bnds)))))
1338                             (,%begin (,cfname . ,(cdr fn))
1339                                      (values . ,ret)))))))
1340)
1341
1342     
1343(define-syntax blas-level3-wrapx 
1344   (lambda (x r c)
1345     (let* ((fn     (cadr x))
1346            (ret    (caddr x))
1347            (errs   (cadddr x)))
1348       
1349       `(begin
1350         (blas-level3-wrap ,(cons (string->symbol (conc "s" (symbol->string (car fn)))) (cdr fn))
1351                           ,ret ,errs #f #f)
1352         (blas-level3-wrap ,(cons (string->symbol (conc "d" (symbol->string (car fn)))) (cdr fn))
1353                           ,ret ,errs #f #f)
1354         (blas-level3-wrap ,(cons (string->symbol (conc "c" (symbol->string (car fn)))) (cdr fn))
1355                           ,ret ,errs #f #f)
1356         (blas-level3-wrap ,(cons (string->symbol (conc "z" (symbol->string (car fn)))) (cdr fn))
1357                           ,ret ,errs #f #f)
1358         
1359         (blas-level3-wrap ,(cons (string->symbol (conc "s" (symbol->string (car fn)))) (cdr fn))
1360                           ,ret ,errs f32vector-length #f)
1361         (blas-level3-wrap ,(cons (string->symbol (conc "d" (symbol->string (car fn)))) (cdr fn))
1362                           ,ret ,errs f64vector-length #f)
1363         (blas-level3-wrap ,(cons (string->symbol (conc "c" (symbol->string (car fn)))) (cdr fn))
1364                           ,ret ,errs (lambda (v) (fx/ 2 (f32vector-length v))) #f)
1365         (blas-level3-wrap ,(cons (string->symbol (conc "z" (symbol->string (car fn)))) (cdr fn))
1366                           ,ret ,errs (lambda (v) (fx/ 2 (f64vector-length v))) #f)
1367         
1368         (blas-level3-wrap ,(cons (string->symbol (conc "s" (symbol->string (car fn)))) (cdr fn))
1369                           ,ret ,errs f32vector-length  scopy)
1370         (blas-level3-wrap ,(cons (string->symbol (conc "d" (symbol->string (car fn)))) (cdr fn))
1371                           ,ret ,errs f64vector-length  dcopy)
1372         (blas-level3-wrap ,(cons (string->symbol (conc "c" (symbol->string (car fn)))) (cdr fn))
1373                           ,ret ,errs (lambda (v) (fx/ 2 (f32vector-length v))) ccopy)
1374         (blas-level3-wrap ,(cons (string->symbol (conc "z" (symbol->string (car fn)))) (cdr fn))
1375                           ,ret ,errs (lambda (v) (fx/ 2 (f64vector-length v))) zcopy))))
1376
1377   )
1378
1379
1380(define-syntax blas-level3-cz-wrapx
1381   (lambda (x r c)
1382     (let* ((fn      (cadr x))
1383            (ret     (caddr x))
1384            (errs    (cadddr x)))
1385       
1386       `(begin
1387         (blas-level3-wrap ,(cons (string->symbol (conc "c" (symbol->string (car fn)))) (cdr fn))
1388                           ,ret ,errs #f #f)
1389         (blas-level3-wrap ,(cons (string->symbol (conc "z" (symbol->string (car fn)))) (cdr fn))
1390                           ,ret ,errs #f #f)
1391         
1392         (blas-level3-wrap ,(cons (string->symbol (conc "c" (symbol->string (car fn)))) (cdr fn))
1393                           ,ret ,errs (lambda (v) (fx/ 2 (f32vector-length v))) #f)
1394         (blas-level3-wrap ,(cons (string->symbol (conc "z" (symbol->string (car fn)))) (cdr fn))
1395                           ,ret ,errs (lambda (v) (fx/ 2 (f64vector-length v))) #f)
1396         
1397         (blas-level3-wrap ,(cons (string->symbol (conc "c" (symbol->string (car fn)))) (cdr fn))
1398                           ,ret ,errs (lambda (v) (fx/ 2 (f32vector-length v))) ccopy)
1399         (blas-level3-wrap ,(cons (string->symbol (conc "z" (symbol->string (car fn)))) (cdr fn))
1400                           ,ret ,errs (lambda (v) (fx/ 2 (f64vector-length v))) zcopy))))
1401  )
1402
1403
1404(blas-level3-wrapx (gemm order transA transB m n k alpha a lda b ldb beta c ldc)
1405                     (c)
1406                     (lambda (i) (cond ((= i 3)  "M < 0")
1407                                       ((= i 4)  "N < 0")
1408                                       ((= i 5)  "K < 0")
1409                                       ((= i 8)  "LDA < max(1, M or K)")
1410                                       ((= i 10) "LDB < max(1, N or K)")
1411                                       ((= i 13) "LDC < max(1, M)")
1412                                       (else (conc "error code " i)))))
1413 
1414  (blas-level3-wrapx (symm order side uplo  m n alpha a lda b ldb beta c ldc)
1415                     (c)
1416                     (lambda (i) (cond ((= i 3)  "M < 0")
1417                                     ((= i 4)  "N < 0")
1418                                     ((= i 5)  "K < 0")
1419                                     ((= i 8)  "LDA < max(1, M or K)")
1420                                     ((= i 10) "LDB < max(1, N or K)")
1421                                     ((= i 13) "LDC < max(1, M)")
1422                                     (else (conc "error code " i)))))
1423 
1424  (blas-level3-cz-wrapx (hemm order side uplo  m n alpha a lda b ldb beta c ldc)
1425                        (c)
1426                        (lambda (i) (cond ((= i 3)  "M < 0")
1427                                          ((= i 4)  "N < 0")
1428                                          ((= i 5)  "K < 0")
1429                                          ((= i 8)  "LDA < max(1, M or K)")
1430                                          ((= i 10) "LDB < max(1, N or K)")
1431                                          ((= i 13) "LDC < max(1, M)")
1432                                          (else (conc "error code " i)))))
1433 
1434  (blas-level3-wrapx (syrk order uplo trans n k alpha a lda beta c ldc)
1435                     (c)
1436                     (lambda (i) (cond ((= i 3)  "M < 0")
1437                                       ((= i 4)  "N < 0")
1438                                       ((= i 5)  "K < 0")
1439                                       ((= i 8)  "LDA < max(1, M or K)")
1440                                       ((= i 10) "LDB < max(1, N or K)")
1441                                       ((= i 13) "LDC < max(1, M)")
1442                                       (else (conc "error code " i)))))
1443 
1444  (blas-level3-cz-wrapx (herk order uplo trans n k alpha a lda beta c ldc)
1445                        (c)
1446                        (lambda (i) (cond ((= i 3)  "M < 0")
1447                                          ((= i 4)  "N < 0")
1448                                          ((= i 5)  "K < 0")
1449                                          ((= i 8)  "LDA < max(1, M or K)")
1450                                          ((= i 10) "LDB < max(1, N or K)")
1451                                          ((= i 13) "LDC < max(1, M)")
1452                                          (else (conc "error code " i)))))
1453 
1454 
1455  (blas-level3-wrapx (syr2k order uplo trans n k alpha a lda b ldb beta c ldc)
1456                     (c)
1457                     (lambda (i) (cond ((= i 3)  "M < 0")
1458                                       ((= i 4)  "N < 0")
1459                                       ((= i 5)  "K < 0")
1460                                       ((= i 8)  "LDA < max(1, M or K)")
1461                                       ((= i 10) "LDB < max(1, N or K)")
1462                                       ((= i 13) "LDC < max(1, M)")
1463                                       (else (conc "error code " i)))))
1464 
1465  (blas-level3-cz-wrapx (her2k order uplo trans n k alpha a lda b ldb beta c ldc)
1466                        (c)
1467                        (lambda (i) (cond ((= i 3)  "M < 0")
1468                                          ((= i 4)  "N < 0")
1469                                          ((= i 5)  "K < 0")
1470                                          ((= i 8)  "LDA < max(1, M or K)")
1471                                          ((= i 10) "LDB < max(1, N or K)")
1472                                          ((= i 13) "LDC < max(1, M)")
1473                                          (else (conc "error code " i)))))
1474 
1475  (blas-level3-wrapx (trmm order side uplo transA diag m n alpha a lda b ldb)
1476                     (b)
1477                     (lambda (i) (cond ((= i 3)  "M < 0")
1478                                       ((= i 4)  "N < 0")
1479                                       ((= i 5)  "K < 0")
1480                                       ((= i 8)  "LDA < max(1, M or K)")
1481                                       ((= i 10) "LDB < max(1, N or K)")
1482                                       ((= i 13) "LDC < max(1, M)")
1483                                       (else (conc "error code " i)))))
1484 
1485 
1486  (blas-level3-wrapx (trsm order side uplo transA diag m n alpha a lda b ldb)
1487                     (b)
1488                     (lambda (i) (cond ((= i 3)  "M < 0")
1489                                       ((= i 4)  "N < 0")
1490                                       ((= i 5)  "K < 0")
1491                                       ((= i 8)  "LDA < max(1, M or K)")
1492                                       ((= i 10) "LDB < max(1, N or K)")
1493                                       ((= i 13) "LDC < max(1, M)")
1494                                       (else (conc "error code " i)))))
1495     
1496
1497
1498(define-syntax blas-level2-wrap
1499      (lambda (x r c)
1500        (let* ((fn      (cadr x))
1501               (ret     (caddr x))
1502               (err     (cadddr x))
1503               (vsize   (car (cddddr x)))
1504               (copy    (cadr (cddddr x)))
1505               (cfname  (string->symbol (conc "cblas_" (symbol->string (car fn)))))
1506               (fname   (string->symbol (conc (if vsize "" "unsafe-") 
1507                                              (symbol->string (car fn))
1508                                              (if copy "" "!"))))
1509               (%define         (r 'define))
1510               (%begin          (r 'begin))
1511               (%let            (r 'let))
1512               (%cond           (r 'cond))
1513               (%or             (r 'or))
1514               (%if             (r 'if))
1515               (%let-optionals  (r 'let-optionals))
1516
1517               (ka              (r 'ka))
1518               (asize           (r 'asize))
1519               (apsize          (r 'apsize))
1520               (apdim           (r 'apdim))
1521               (xsize           (r 'xsize))
1522               (ysize           (r 'ysize))
1523               (xdim            (r 'xdim))
1524               (ydim            (r 'ydim))
1525
1526               (args  (reverse (cdr fn)))
1527
1528               (fsig  (let loop ((args args) (sig 'rest))
1529                        (if (null? args) (cons fname sig)
1530                            (let ((x (car args)))
1531                              (let ((sig (case x
1532                                          ((lda)      sig)
1533                                          ((incx)     sig)
1534                                          ((incy)     sig)
1535                                          ((offx)     sig)
1536                                          ((offy)     sig)
1537                                          (else      (cons x sig)))))
1538                                (loop (cdr args) sig))))))
1539
1540               (opts  (append
1541                       (if (memq 'lda fn)  `((lda  ,(cond ((memq 'k fn) `(fx+ 1 k))
1542                                                          (else 'n)))) `())
1543                       (if (memq 'incy fn) `((incx 1) (incy 1) (offx 0) (offy 0)) `((incx 1)))))
1544               )
1545
1546      `(,%define ,fsig 
1547         (,%let-optionals rest ,opts
1548          ,(if vsize
1549               `(,%begin
1550                  ,(if (memq 'a fn)
1551                       `(,%let ((,asize (,vsize a))
1552                                (,ka    ,(if (memq 'm fn) 'm 'n)))
1553                          (,%if (< ,asize (fx* lda ,ka))
1554                              (blas:error ',fname (conc "matrix A is allocated " ,asize " elements "
1555                                                        "but given dimensions are " ,ka " by " lda))))
1556                       `(begin))
1557                  ,(if (memq 'ap fn)
1558                       `(,%let ((,apsize (,vsize ap))
1559                                (,apdim  (fx/ (fx* n (fx+ n 1)) 2)))
1560                          (,%if (< ,apsize ,apdim)
1561                                (blas:error ',fname (conc "vector Ap is allocated " ,apsize " elements "
1562                                                          "but given dimension is " ,apdim))))
1563                       `(begin))
1564                  ,(if (memq 'y fn)
1565                       `(,%let ((,ysize (,vsize y))
1566                                (,ydim  ,(if (and (memq 'm fn) (memq 'trans fn))
1567                                             `(,%if (= trans NoTrans) 
1568                                                    (fx+ 1 (fx* (abs incy) (fx- (fx+ offy m) 1)))
1569                                                    (fx+ 1 (fx* (abs incy) (fx- (fx+ offy n) 1))))
1570                                             `(fx+ 1 (fx* (abs incy) (fx- n 1))))))
1571                          (,%if (< ,ysize ,ydim)
1572                                (blas:error ',fname (conc "vector Y is allocated " ,ysize " elements "
1573                                                          "but given dimension is " ,ydim))))
1574                       `(begin))
1575                  ,(if (memq 'x fn)
1576                       `(,%let ((,xsize (,vsize x))
1577                                (,xdim  ,(if (and (memq 'm fn) (memq 'trans fn))
1578                                             `(if (= trans NoTrans) 
1579                                                  (fx+ 1 (fx* (abs incx) (fx- (fx+ offx n) 1)))
1580                                                  (fx+ 1 (fx* (abs incx) (fx- (fx+ offx m) 1))))
1581                                             `(fx+ 1 (fx* (abs incx) (fx- n 1))))))
1582                          (,%if (< ,xsize ,xdim)
1583                                (blas:error ',fname (conc "vector X is allocated " ,xsize " elements "
1584                                                          "but given dimension is " ,xdim))))
1585                       `(begin)))
1586               `(begin))
1587          (let ,(let loop ((fn fn) (bnds '()))
1588                  (if (null? fn) bnds
1589                      (let ((x (car fn)))
1590                        (let ((bnds (case x 
1591                                      (else    (if (and copy (memq x ret))
1592                                                   (cons `(,x (,copy ,x)) bnds)
1593                                                   bnds)))))
1594                          (loop (cdr fn) bnds)))))
1595            (begin (,cfname . ,(cdr fn))
1596                   (values . ,ret)))))))
1597)
1598
1599
1600(define-syntax blas-level2-wrapx
1601      (lambda (x r c)
1602        (let* ((fn      (cadr x))
1603               (ret     (caddr x))
1604               (errs    (cadddr x)))
1605         
1606          `(begin
1607            (blas-level2-wrap ,(cons (string->symbol (conc "s" (symbol->string (car fn)))) (cdr fn))
1608                              ,ret ,errs #f #f)
1609            (blas-level2-wrap ,(cons (string->symbol (conc "d" (symbol->string (car fn)))) (cdr fn))
1610                              ,ret ,errs #f #f)
1611            (blas-level2-wrap ,(cons (string->symbol (conc "c" (symbol->string (car fn)))) (cdr fn))
1612                              ,ret ,errs #f #f)
1613            (blas-level2-wrap ,(cons (string->symbol (conc "z" (symbol->string (car fn)))) (cdr fn))
1614                              ,ret ,errs #f #f)
1615           
1616            (blas-level2-wrap ,(cons (string->symbol (conc "s" (symbol->string (car fn)))) (cdr fn))
1617                              ,ret ,errs f32vector-length #f)
1618            (blas-level2-wrap ,(cons (string->symbol (conc "d" (symbol->string (car fn)))) (cdr fn))
1619                              ,ret ,errs f64vector-length #f)
1620            (blas-level2-wrap ,(cons (string->symbol (conc "c" (symbol->string (car fn)))) (cdr fn))
1621                              ,ret ,errs (lambda (v) (fx/ 2 (f32vector-length v))) #f)
1622            (blas-level2-wrap ,(cons (string->symbol (conc "z" (symbol->string (car fn)))) (cdr fn))
1623                              ,ret ,errs (lambda (v) (fx/ 2 (f64vector-length v))) #f)
1624           
1625            (blas-level2-wrap ,(cons (string->symbol (conc "s" (symbol->string (car fn)))) (cdr fn))
1626                              ,ret ,errs f32vector-length  scopy)
1627            (blas-level2-wrap ,(cons (string->symbol (conc "d" (symbol->string (car fn)))) (cdr fn))
1628                              ,ret ,errs f64vector-length  dcopy)
1629            (blas-level2-wrap ,(cons (string->symbol (conc "c" (symbol->string (car fn)))) (cdr fn))
1630                              ,ret ,errs (lambda (v) (fx/ 2 (f32vector-length v))) ccopy)
1631            (blas-level2-wrap ,(cons (string->symbol (conc "z" (symbol->string (car fn)))) (cdr fn))
1632                              ,ret ,errs (lambda (v) (fx/ 2 (f64vector-length v))) zcopy)))
1633        ))
1634
1635     
1636(define-syntax blas-level2-sd-wrapx
1637      (lambda (x r c)
1638        (let* ((fn      (cadr x))
1639               (ret     (caddr x))
1640               (errs    (cadddr x)))
1641         
1642          `(begin
1643            (blas-level2-wrap ,(cons (string->symbol (conc "s" (symbol->string (car fn)))) (cdr fn))
1644                              ,ret ,errs #f #f)
1645            (blas-level2-wrap ,(cons (string->symbol (conc "d" (symbol->string (car fn)))) (cdr fn))
1646                              ,ret ,errs #f #f)
1647           
1648            (blas-level2-wrap ,(cons (string->symbol (conc "s" (symbol->string (car fn)))) (cdr fn))
1649                              ,ret ,errs f32vector-length #f)
1650            (blas-level2-wrap ,(cons (string->symbol (conc "d" (symbol->string (car fn)))) (cdr fn))
1651                              ,ret ,errs f64vector-length #f)
1652           
1653            (blas-level2-wrap ,(cons (string->symbol (conc "s" (symbol->string (car fn)))) (cdr fn))
1654                              ,ret ,errs f32vector-length  scopy)
1655            (blas-level2-wrap ,(cons (string->symbol (conc "d" (symbol->string (car fn)))) (cdr fn))
1656                              ,ret ,errs f64vector-length  dcopy))))
1657)
1658 
1659   
1660(define-syntax blas-level2-cz-wrapx
1661      (lambda (x r c)
1662        (let* ((fn      (cadr x))
1663               (ret     (caddr x))
1664               (errs    (cadddr x)))
1665         
1666          `(begin
1667            (blas-level2-wrap ,(cons (string->symbol (conc "c" (symbol->string (car fn)))) (cdr fn))
1668                              ,ret ,errs #f #f)
1669            (blas-level2-wrap ,(cons (string->symbol (conc "z" (symbol->string (car fn)))) (cdr fn))
1670                              ,ret ,errs #f #f)
1671           
1672            (blas-level2-wrap ,(cons (string->symbol (conc "c" (symbol->string (car fn)))) (cdr fn))
1673                              ,ret ,errs (lambda (v) (fx/ 2 (f32vector-length v))) #f)
1674            (blas-level2-wrap ,(cons (string->symbol (conc "z" (symbol->string (car fn)))) (cdr fn))
1675                              ,ret ,errs (lambda (v) (fx/ 2 (f64vector-length v))) #f)
1676           
1677            (blas-level2-wrap ,(cons (string->symbol (conc "c" (symbol->string (car fn)))) (cdr fn))
1678                              ,ret ,errs (lambda (v) (fx/ 2 (f32vector-length v))) ccopy)
1679            (blas-level2-wrap ,(cons (string->symbol (conc "z" (symbol->string (car fn)))) (cdr fn))
1680                              ,ret ,errs (lambda (v) (fx/ 2 (f64vector-length v))) zcopy))))
1681)
1682     
1683     
1684(blas-level2-wrapx (gemv order trans m n alpha a lda x incx beta y incy)
1685                     (y)
1686                     (lambda (i) (cond ((= i 2)  "M < 0")
1687                                       ((= i 3)  "N < 0")
1688                                       ((= i 6)  "LDA < max(1, M)")
1689                                       ((= i 8)  "INCX = 0")
1690                                       ((= i 11) "INCY < = 0")
1691                                       (else (conc "error code " i)))))
1692 
1693(blas-level2-cz-wrapx (hemv order uplo n alpha a lda x incx beta y incy)
1694                      (y)
1695                      (lambda (i) (cond ((= i 2)  "M < 0")
1696                                        ((= i 3)  "N < 0")
1697                                        ((= i 6)  "LDA < max(1, M)")
1698                                        ((= i 8)  "INCX = 0")
1699                                        ((= i 11) "INCY < = 0")
1700                                        (else (conc "error code " i)))))
1701
1702(blas-level2-cz-wrapx (hbmv order uplo n k alpha a lda x incx beta y incy)
1703                      (y)
1704                      (lambda (i) (cond ((= i 2)  "M < 0")
1705                                        ((= i 3)  "N < 0")
1706                                        ((= i 6)  "LDA < max(1, M)")
1707                                        ((= i 8)  "INCX = 0")
1708                                        ((= i 11) "INCY < = 0")
1709                                        (else (conc "error code " i)))))
1710
1711(blas-level2-cz-wrapx (hpmv order uplo n alpha ap x incx beta y incy)
1712                      (y)
1713                      (lambda (i) (cond ((= i 2)  "M < 0")
1714                                        ((= i 3)  "N < 0")
1715                                        ((= i 6)  "LDA < max(1, M)")
1716                                        ((= i 8)  "INCX = 0")
1717                                        ((= i 11) "INCY < = 0")
1718                                        (else (conc "error code " i)))))
1719
1720(blas-level2-sd-wrapx (symv order uplo n alpha a lda x incx beta y incy)
1721                      (y)
1722                      (lambda (i) (cond ((= i 2)  "M < 0")
1723                                        ((= i 3)  "N < 0")
1724                                        ((= i 6)  "LDA < max(1, M)")
1725                                        ((= i 8)  "INCX = 0")
1726                                        ((= i 11) "INCY < = 0")
1727                                        (else (conc "error code " i)))))
1728
1729(blas-level2-sd-wrapx (sbmv order uplo n k alpha a lda x incx beta y incy)
1730                      (y)
1731                      (lambda (i) (cond ((= i 2)  "M < 0")
1732                                        ((= i 3)  "N < 0")
1733                                        ((= i 6)  "LDA < max(1, M)")
1734                                        ((= i 8)  "INCX = 0")
1735                                        ((= i 11) "INCY < = 0")
1736                                        (else (conc "error code " i)))))
1737
1738(blas-level2-sd-wrapx (spmv order uplo n alpha ap x incx beta y incy)
1739                      (y)
1740                      (lambda (i) (cond ((= i 2)  "M < 0")
1741                                        ((= i 3)  "N < 0")
1742                                        ((= i 6)  "LDA < max(1, M)")
1743                                        ((= i 8)  "INCX = 0")
1744                                        ((= i 11) "INCY < = 0")
1745                                        (else (conc "error code " i)))))
1746
1747(blas-level2-wrapx (trmv order uplo trans diag n a lda x incx)
1748                   (x)
1749                   (lambda (i) (cond ((= i 2)  "M < 0")
1750                                     ((= i 3)  "N < 0")
1751                                     ((= i 6)  "LDA < max(1, M)")
1752                                     ((= i 8)  "INCX = 0")
1753                                     ((= i 11) "INCY < = 0")
1754                                     (else (conc "error code " i)))))
1755
1756(blas-level2-wrapx (tbmv order uplo trans diag n k a lda x incx)
1757                   (x)
1758                   (lambda (i) (cond ((= i 2)  "M < 0")
1759                                     ((= i 3)  "N < 0")
1760                                     ((= i 6)  "LDA < max(1, M)")
1761                                     ((= i 8)  "INCX = 0")
1762                                     ((= i 11) "INCY < = 0")
1763                                     (else (conc "error code " i)))))
1764
1765(blas-level2-wrapx (tpmv order uplo trans diag n ap x incx)
1766                   (x)
1767                   (lambda (i) (cond ((= i 2)  "M < 0")
1768                                     ((= i 3)  "N < 0")
1769                                     ((= i 6)  "LDA < max(1, M)")
1770                                     ((= i 8)  "INCX = 0")
1771                                     ((= i 11) "INCY < = 0")
1772                                     (else (conc "error code " i)))))
1773
1774(blas-level2-wrapx (trsv order uplo trans diag n a lda x incx)
1775                   (x)
1776                   (lambda (i) (cond ((= i 2)  "M < 0")
1777                                     ((= i 3)  "N < 0")
1778                                     ((= i 6)  "LDA < max(1, M)")
1779                                     ((= i 8)  "INCX = 0")
1780                                     ((= i 11) "INCY < = 0")
1781                                     (else (conc "error code " i)))))
1782
1783  (blas-level2-wrapx (tbsv order uplo trans diag n k a lda x incx)
1784                     (x)
1785                     (lambda (i) (cond ((= i 2)  "M < 0")
1786                                       ((= i 3)  "N < 0")
1787                                       ((= i 6)  "LDA < max(1, M)")
1788                                       ((= i 8)  "INCX = 0")
1789                                       ((= i 11) "INCY < = 0")
1790                                       (else (conc "error code " i)))))
1791 
1792  (blas-level2-wrapx (tpsv order uplo trans diag n ap x incx)
1793                     (x)
1794                     (lambda (i) (cond ((= i 2)  "M < 0")
1795                                       ((= i 3)  "N < 0")
1796                                       ((= i 6)  "LDA < max(1, M)")
1797                                       ((= i 8)  "INCX = 0")
1798                                       ((= i 11) "INCY < = 0")
1799                                       (else (conc "error code " i)))))
1800 
1801  (blas-level2-sd-wrapx (ger order m n alpha x incx y incy a lda)
1802                        (a)
1803                        (lambda (i) (cond ((= i 2)  "M < 0")
1804                                          ((= i 3)  "N < 0")
1805                                          ((= i 6)  "LDA < max(1, M)")
1806                                          ((= i 8)  "INCX = 0")
1807                                          ((= i 11) "INCY < = 0")
1808                                          (else (conc "error code " i)))))
1809 
1810  (blas-level2-cz-wrapx (geru order m n alpha x incx y incy a lda)
1811                        (a)
1812                        (lambda (i) (cond ((= i 2)  "M < 0")
1813                                          ((= i 3)  "N < 0")
1814                                          ((= i 6)  "LDA < max(1, M)")
1815                                          ((= i 8)  "INCX = 0")
1816                                          ((= i 11) "INCY < = 0")
1817                                          (else (conc "error code " i)))))
1818                             
1819  (blas-level2-cz-wrapx (gerc order m n alpha x incx y incy a lda)
1820                        (a)
1821                        (lambda (i) (cond ((= i 2)  "M < 0")
1822                                          ((= i 3)  "N < 0")
1823                                          ((= i 6)  "LDA < max(1, M)")
1824                                          ((= i 8)  "INCX = 0")
1825                                          ((= i 11) "INCY < = 0")
1826                                          (else (conc "error code " i)))))
1827 
1828                             
1829  (blas-level2-cz-wrapx (her order uplo n alpha x incx a lda)
1830                        (a)
1831                        (lambda (i) (cond ((= i 2)  "M < 0")
1832                                          ((= i 3)  "N < 0")
1833                                          ((= i 6)  "LDA < max(1, M)")
1834                                          ((= i 8)  "INCX = 0")
1835                                          ((= i 11) "INCY < = 0")
1836                                          (else (conc "error code " i)))))
1837 
1838  (blas-level2-cz-wrapx (hpr order uplo n alpha x incx ap)
1839                        (ap)
1840                        (lambda (i) (cond ((= i 2)  "M < 0")
1841                                          ((= i 3)  "N < 0")
1842                                          ((= i 6)  "LDA < max(1, M)")
1843                                          ((= i 8)  "INCX = 0")
1844                                          ((= i 11) "INCY < = 0")
1845                                          (else (conc "error code " i)))))
1846 
1847  (blas-level2-cz-wrapx (her2 order uplo n alpha x incx y incy a lda)
1848                        (a)
1849                        (lambda (i) (cond ((= i 2)  "M < 0")
1850                                          ((= i 3)  "N < 0")
1851                                          ((= i 6)  "LDA < max(1, M)")
1852                                          ((= i 8)  "INCX = 0")
1853                                          ((= i 11) "INCY < = 0")
1854                                          (else (conc "error code " i)))))
1855 
1856  (blas-level2-cz-wrapx (hpr2 order uplo n alpha x incx y incy ap)
1857                        (ap)
1858                        (lambda (i) (cond ((= i 2)  "M < 0")
1859                                          ((= i 3)  "N < 0")
1860                                          ((= i 6)  "LDA < max(1, M)")
1861                                          ((= i 8)  "INCX = 0")
1862                                          ((= i 11) "INCY < = 0")
1863                                          (else (conc "error code " i)))))
1864 
1865  (blas-level2-sd-wrapx (syr order uplo n alpha x incx a lda)
1866                        (a)
1867                        (lambda (i) (cond ((= i 2)  "M < 0")
1868                                          ((= i 3)  "N < 0")
1869                                          ((= i 6)  "LDA < max(1, M)")
1870                                          ((= i 8)  "INCX = 0")
1871                                          ((= i 11) "INCY < = 0")
1872                                          (else (conc "error code " i)))))
1873                             
1874  (blas-level2-sd-wrapx (spr order uplo n alpha x incx ap)
1875                        (ap)
1876                        (lambda (i) (cond ((= i 2)  "M < 0")
1877                                          ((= i 3)  "N < 0")
1878                                          ((= i 6)  "LDA < max(1, M)")
1879                                          ((= i 8)  "INCX = 0")
1880                                          ((= i 11) "INCY < = 0")
1881                                          (else (conc "error code " i)))))
1882 
1883  (blas-level2-sd-wrapx (syr2 order uplo n alpha x incx y incy a lda)
1884                        (a)
1885                        (lambda (i) (cond ((= i 2)  "M < 0")
1886                                          ((= i 3)  "N < 0")
1887                                          ((= i 6)  "LDA < max(1, M)")
1888                                          ((= i 8)  "INCX = 0")
1889                                          ((= i 11) "INCY < = 0")
1890                                          (else (conc "error code " i)))))
1891 
1892  (blas-level2-sd-wrapx (ger order m n alpha x incx y incy a lda)
1893                        (a)
1894                        (lambda (i) (cond ((= i 2)  "M < 0")
1895                                          ((= i 3)  "N < 0")
1896                                          ((= i 6)  "LDA < max(1, M)")
1897                                          ((= i 8)  "INCX = 0")
1898                                          ((= i 11) "INCY < = 0")
1899                                          (else (conc "error code " i)))))
1900 
1901  (blas-level2-sd-wrapx (iger order m n alpha x incx offx y incy offy a lda)
1902                        (a)
1903                        (lambda (i) (cond ((= i 2)  "M < 0")
1904                                          ((= i 3)  "N < 0")
1905                                          ((= i 6)  "LDA < max(1, M)")
1906                                          ((= i 8)  "INCX = 0")
1907                                          ((= i 11) "INCY < = 0")
1908                                          (else (conc "error code " i)))))
1909 
1910  (blas-level2-cz-wrapx (geru order m n alpha x incx y incy a lda)
1911                        (a)
1912                        (lambda (i) (cond ((= i 2)  "M < 0")
1913                                          ((= i 3)  "N < 0")
1914                                          ((= i 6)  "LDA < max(1, M)")
1915                                          ((= i 8)  "INCX = 0")
1916                                          ((= i 11) "INCY < = 0")
1917                                          (else (conc "error code " i)))))
1918 
1919  (blas-level2-cz-wrapx (gerc order m n alpha x incx y incy a lda)
1920                        (a)
1921                        (lambda (i) (cond ((= i 2)  "M < 0")
1922                                          ((= i 3)  "N < 0")
1923                                          ((= i 6)  "LDA < max(1, M)")
1924                                          ((= i 8)  "INCX = 0")
1925                                          ((= i 11) "INCY < = 0")
1926                                          (else (conc "error code " i)))))
1927 
1928 
1929  (blas-level2-cz-wrapx (her order uplo n alpha x incx a lda)
1930                        (a)
1931                        (lambda (i) (cond ((= i 2)  "M < 0")
1932                                          ((= i 3)  "N < 0")
1933                                          ((= i 6)  "LDA < max(1, M)")
1934                                          ((= i 8)  "INCX = 0")
1935                                          ((= i 11) "INCY < = 0")
1936                                          (else (conc "error code " i)))))
1937 
1938  (blas-level2-cz-wrapx (hpr order uplo n alpha x incx ap)
1939                        (ap)
1940                        (lambda (i) (cond ((= i 2)  "M < 0")
1941                                          ((= i 3)  "N < 0")
1942                                          ((= i 6)  "LDA < max(1, M)")
1943                                          ((= i 8)  "INCX = 0")
1944                                          ((= i 11) "INCY < = 0")
1945                                          (else (conc "error code " i)))))
1946 
1947  (blas-level2-cz-wrapx (her2 order uplo n alpha x incx y incy a lda)
1948                        (a)
1949                        (lambda (i) (cond ((= i 2)  "M < 0")
1950                                          ((= i 3)  "N < 0")
1951                                          ((= i 6)  "LDA < max(1, M)")
1952                                          ((= i 8)  "INCX = 0")
1953                                          ((= i 11) "INCY < = 0")
1954                                          (else (conc "error code " i)))))
1955 
1956  (blas-level2-cz-wrapx (hpr2 order uplo n alpha x incx y incy ap)
1957                        (ap)
1958                        (lambda (i) (cond ((= i 2)  "M < 0")
1959                                          ((= i 3)  "N < 0")
1960                                          ((= i 6)  "LDA < max(1, M)")
1961                                          ((= i 8)  "INCX = 0")
1962                                          ((= i 11) "INCY < = 0")
1963                                          (else (conc "error code " i)))))
1964                             
1965  (blas-level2-sd-wrapx (syr order uplo n alpha x incx a lda)
1966                        (a)
1967                        (lambda (i) (cond ((= i 2)  "M < 0")
1968                                          ((= i 3)  "N < 0")
1969                                          ((= i 6)  "LDA < max(1, M)")
1970                                          ((= i 8)  "INCX = 0")
1971                                          ((= i 11) "INCY < = 0")
1972                                          (else (conc "error code " i)))))
1973 
1974  (blas-level2-sd-wrapx (spr order uplo n alpha x incx ap)
1975                        (ap)
1976                        (lambda (i) (cond ((= i 2)  "M < 0")
1977                                          ((= i 3)  "N < 0")
1978                                          ((= i 6)  "LDA < max(1, M)")
1979                                          ((= i 8)  "INCX = 0")
1980                                          ((= i 11) "INCY < = 0")
1981                                          (else (conc "error code " i)))))
1982 
1983  (blas-level2-sd-wrapx (syr2 order uplo n alpha x incx y incy a lda)
1984                        (a)
1985                        (lambda (i) (cond ((= i 2)  "M < 0")
1986                                          ((= i 3)  "N < 0")
1987                                          ((= i 6)  "LDA < max(1, M)")
1988                                          ((= i 8)  "INCX = 0")
1989                                          ((= i 11) "INCY < = 0")
1990                                          (else (conc "error code " i)))))
1991 
1992  (blas-level2-sd-wrapx (spr2 order uplo n alpha x incx y incy ap)
1993                        (ap)
1994                        (lambda (i) (cond ((= i 2)  "M < 0")
1995                                          ((= i 3)  "N < 0")
1996                                          ((= i 6)  "LDA < max(1, M)")
1997                                          ((= i 8)  "INCX = 0")
1998                                          ((= i 11) "INCY < = 0")
1999                                          (else (conc "error code " i)))))
2000 
2001   
2002
2003
2004(define-syntax  blas-level1-wrap
2005      (lambda (x r c)
2006        (let* ((fn            (cadr x))
2007               (ret           (caddr x))
2008               (err           (cadddr x))
2009               (vsize         (car (cddddr x)))
2010               (copy          (cadr (cddddr x)))
2011               (make-return   (cddr (cddddr x)))
2012               (cfname  (string->symbol (conc "cblas_" (symbol->string (car fn)))))
2013               (fname   (string->symbol (conc (if vsize "" "unsafe-") 
2014                                              (symbol->string (car fn))
2015                                              (if copy "" "!"))))
2016               
2017               (%define         (r 'define))
2018               (%begin          (r 'begin))
2019               (%let            (r 'let))
2020               (%cond           (r 'cond))
2021               (%or             (r 'or))
2022               (%if             (r 'if))
2023               (%let-optionals  (r 'let-optionals))
2024
2025               (asize           (r 'asize))
2026               (apsize          (r 'apsize))
2027               (apdim           (r 'apdim))
2028               (xsize           (r 'xsize))
2029               (ysize           (r 'ysize))
2030               (xdim            (r 'xdim))
2031               (ydim            (r 'ydim))
2032               (psize           (r 'psize))
2033               (pdim            (r 'pdim))
2034
2035               (args   (reverse (cdr fn)))
2036
2037               (fsig  (let loop ((args args) (sig 'rest))
2038                        (if (null? args) (cons fname sig)
2039                            (let ((x (car args)))
2040                              (let ((sig (case x
2041                                          ((incx)     sig)
2042                                          ((incy)     sig)
2043                                          ((dotu)     sig)
2044                                          ((dotc)     sig)
2045                                          ((offx)     sig)
2046                                          ((offy)     sig)
2047                                          (else      (cons x sig)))))
2048                                (loop (cdr args) sig))))))
2049
2050               (opts  (cond ((memq 'incy fn)  `((incx 1) (incy 1) (offx 0) (offy 0)))
2051                            (else `((incx 1) (offx 0))))))
2052
2053          `(,%define ,fsig 
2054                     (,%let-optionals rest ,opts
2055                      ,(if vsize
2056                           `(,%begin
2057                             ,(if (memq 'y fn)
2058                                  `(,%let ((,ysize (,vsize y))
2059                                           (,ydim  (fx+ 1 (fx* (abs incy) (fx- (fx+ offy n) 1)))))
2060                                          (,%if (< ,ysize ,ydim)
2061                                                (blas:error ',fname (conc "vector Y is allocated " ,ysize " elements "
2062                                                                          "but given dimension is " ,ydim))))
2063                                  `(begin))
2064                             ,(if (memq 'x fn)
2065                                  `(,%let ((,xsize (,vsize x))
2066                                           (,xdim  (fx+ 1 (fx* (abs incx) (fx- (fx+ offx n) 1)))))
2067                                          (,%if (< ,xsize ,xdim)
2068                                                (blas:error ',fname (conc "vector X is allocated " ,xsize " elements "
2069                                                                          "but given dimension is " ,xdim))))
2070                                  `(begin))
2071                             ,(if (memq 'param fn)
2072                                  `(,%let ((,psize (,vsize param))
2073                                           (,pdim  5))
2074                                          (,%if (< ,psize ,pdim)
2075                                                (blas:error ',fname (conc "vector PARAM is allocated " ,psize " elements "
2076                                                                          "but dimension must be " ,pdim))))
2077                                  `(begin)))
2078                           
2079                           `(begin))
2080                      (let ,(let loop ((fn fn) (bnds '()))
2081                              (if (null? fn) bnds
2082                                  (let ((x (car fn)))
2083                                    (let ((bnds (cond ((or (eq? x 'dotc) (eq? x 'dotu))
2084                                                       (cons `(,x (,(car make-return))) bnds))
2085                                                      ((and copy (memq x ret))
2086                                                       (cons `(,x (,copy ,x)) bnds))
2087                                                      (else bnds))))
2088                                      (loop (cdr fn) bnds)))))
2089                        ,(cond
2090                          ((memq 'dotc fn)   `(begin (,cfname . ,(cdr fn))
2091                                                     (values dotc)))
2092                          ((memq 'dotu fn)   `(begin (,cfname . ,(cdr fn))
2093                                                     (values dotu)))
2094                          ((not ret)         `(,cfname . ,(cdr fn)))
2095                          (else              `(begin (,cfname . ,(cdr fn))
2096                                                     (values . ,ret)))))))))
2097      )
2098
2099
2100(define-syntax blas-level1-wrapx
2101       (lambda (x r c) 
2102         (let* ((fn      (cadr x))
2103                (ret     (caddr x))
2104                (errs    (cadddr x)))
2105           
2106           (if (not ret)
2107               `(begin
2108                 (blas-level1-wrap ,(cons (string->symbol (conc "s" (symbol->string (car fn)))) (cdr fn))
2109                               ,ret ,errs f32vector-length  scopy)
2110                 (blas-level1-wrap ,(cons (string->symbol (conc "d" (symbol->string (car fn)))) (cdr fn))
2111                                   ,ret ,errs f64vector-length  dcopy)
2112                 (blas-level1-wrap ,(cons (string->symbol (conc "c" (symbol->string (car fn)))) (cdr fn))
2113                                   ,ret ,errs (lambda (v) (fx/ 2 (f32vector-length v))) ccopy)
2114                 (blas-level1-wrap ,(cons (string->symbol (conc "z" (symbol->string (car fn)))) (cdr fn))
2115                                   ,ret ,errs (lambda (v) (fx/ 2 (f64vector-length v))) zcopy))
2116               
2117               `(begin
2118                 (blas-level1-wrap ,(cons (string->symbol (conc "s" (symbol->string (car fn)))) (cdr fn))
2119                                   ,ret ,errs #f #f)
2120                 (blas-level1-wrap ,(cons (string->symbol (conc "d" (symbol->string (car fn)))) (cdr fn))
2121                                   ,ret ,errs #f #f)
2122                 (blas-level1-wrap ,(cons (string->symbol (conc "c" (symbol->string (car fn)))) (cdr fn))
2123                                   ,ret ,errs #f #f)
2124                 (blas-level1-wrap ,(cons (string->symbol (conc "z" (symbol->string (car fn)))) (cdr fn))
2125                                   ,ret ,errs #f #f)
2126                 
2127                 (blas-level1-wrap ,(cons (string->symbol (conc "s" (symbol->string (car fn)))) (cdr fn))
2128                                   ,ret ,errs f32vector-length #f)
2129                 (blas-level1-wrap ,(cons (string->symbol (conc "d" (symbol->string (car fn)))) (cdr fn))
2130                                   ,ret ,errs f64vector-length #f)
2131                 (blas-level1-wrap ,(cons (string->symbol (conc "c" (symbol->string (car fn)))) (cdr fn))
2132                                   ,ret ,errs (lambda (v) (fx/ 2 (f32vector-length v))) #f)
2133                 (blas-level1-wrap ,(cons (string->symbol (conc "z" (symbol->string (car fn)))) (cdr fn))
2134                                   ,ret ,errs (lambda (v) (fx/ 2 (f64vector-length v))) #f)
2135                 
2136                 (blas-level1-wrap ,(cons (string->symbol (conc "s" (symbol->string (car fn)))) (cdr fn))
2137                                   ,ret ,errs f32vector-length  scopy)
2138                 (blas-level1-wrap ,(cons (string->symbol (conc "d" (symbol->string (car fn)))) (cdr fn))
2139                                   ,ret ,errs f64vector-length  dcopy)
2140                 (blas-level1-wrap ,(cons (string->symbol (conc "c" (symbol->string (car fn)))) (cdr fn))
2141                                   ,ret ,errs (lambda (v) (fx/ 2 (f32vector-length v))) ccopy)
2142                 (blas-level1-wrap ,(cons (string->symbol (conc "z" (symbol->string (car fn)))) (cdr fn))
2143                                   ,ret ,errs (lambda (v) (fx/ 2 (f64vector-length v))) zcopy))))
2144         ))
2145
2146(define-syntax blas-level1-sd-wrapx
2147      (lambda (x r c) 
2148        (let* ((fn      (cadr x))
2149               (ret     (caddr x))
2150               (errs    (cadddr x)))
2151          (if (not ret)
2152             
2153              `(begin
2154                (blas-level1-wrap ,(cons (string->symbol (conc "s" (symbol->string (car fn)))) (cdr fn))
2155                                  ,ret ,errs f32vector-length  scopy)
2156                (blas-level1-wrap ,(cons (string->symbol (conc "d" (symbol->string (car fn)))) (cdr fn))
2157                                  ,ret ,errs f64vector-length  dcopy))
2158             
2159              `(begin
2160                (blas-level1-wrap ,(cons (string->symbol (conc "s" (symbol->string (car fn)))) (cdr fn))
2161                                  ,ret ,errs #f #f)
2162                (blas-level1-wrap ,(cons (string->symbol (conc "d" (symbol->string (car fn)))) (cdr fn))
2163                                  ,ret ,errs #f #f)
2164             
2165                (blas-level1-wrap ,(cons (string->symbol (conc "s" (symbol->string (car fn)))) (cdr fn))
2166                                  ,ret ,errs f32vector-length #f)
2167                (blas-level1-wrap ,(cons (string->symbol (conc "d" (symbol->string (car fn)))) (cdr fn))
2168                                  ,ret ,errs f64vector-length #f)
2169               
2170                (blas-level1-wrap ,(cons (string->symbol (conc "s" (symbol->string (car fn)))) (cdr fn))
2171                                  ,ret ,errs f32vector-length  scopy)
2172                (blas-level1-wrap ,(cons (string->symbol (conc "d" (symbol->string (car fn)))) (cdr fn))
2173                                  ,ret ,errs f64vector-length  dcopy))))
2174        ))
2175
2176
2177(define-syntax blas-level1-cz-wrapx
2178      (lambda (x r c) 
2179        (let* ((fn      (cadr x))
2180               (ret     (caddr x))
2181               (errs    (cadddr x)))
2182         
2183          (if (not ret)
2184              `(begin
2185                (blas-level1-wrap ,(cons (string->symbol (conc "c" (symbol->string (car fn)))) (cdr fn))
2186                                  ,ret ,errs (lambda (v) (fx/ 2 (f32vector-length v))) 
2187                                  ccopy (lambda () (make-f32vector 2)))
2188                (blas-level1-wrap ,(cons (string->symbol (conc "z" (symbol->string (car fn)))) (cdr fn))
2189                                  ,ret ,errs (lambda (v) (fx/ 2 (f64vector-length v))) 
2190                                  zcopy (lambda () (make-f64vector 2))))
2191             
2192              `(begin
2193                (blas-level1-wrap ,(cons (string->symbol (conc "c" (symbol->string (car fn)))) (cdr fn))
2194                                  ,ret ,errs #f #f)
2195                (blas-level1-wrap ,(cons (string->symbol (conc "z" (symbol->string (car fn)))) (cdr fn))
2196                                  ,ret ,errs #f #f)
2197               
2198                (blas-level1-wrap ,(cons (string->symbol (conc "c" (symbol->string (car fn)))) (cdr fn))
2199                                  ,ret ,errs (lambda (v) (fx/ 2 (f32vector-length v))) #f)
2200                (blas-level1-wrap ,(cons (string->symbol (conc "z" (symbol->string (car fn)))) (cdr fn))
2201                                  ,ret ,errs (lambda (v) (fx/ 2 (f64vector-length v))) #f)
2202               
2203                (blas-level1-wrap ,(cons (string->symbol (conc "c" (symbol->string (car fn)))) (cdr fn))
2204                                  ,ret ,errs (lambda (v) (fx/ 2 (f32vector-length v))) ccopy)
2205                (blas-level1-wrap ,(cons (string->symbol (conc "z" (symbol->string (car fn)))) (cdr fn))
2206                                  ,ret ,errs (lambda (v) (fx/ 2 (f64vector-length v))) zcopy))))
2207        ))
2208     
2209     
2210  (blas-level1-sd-wrapx (rot n x incx y incy c s)
2211                        (x y)
2212                        (lambda (i) (cond (conc "error code " i))))
2213 
2214  (blas-level1-sd-wrapx (rotm n x incx y incy param)
2215                        (x y)
2216                        (lambda (i) (cond (conc "error code " i))))
2217 
2218  (blas-level1-wrapx (swap n x incx y incy)
2219                     (x y)
2220                     (lambda (i) (cond (conc "error code " i))))
2221 
2222  (blas-level1-wrapx (scal n alpha x incx)
2223                     (x)
2224                     (lambda (i) (cond (conc "error code " i))))
2225 
2226  (blas-level1-wrapx (axpy n alpha x incx y incy)
2227                     (y)
2228                     (lambda (i) (cond (conc "error code " i))))
2229 
2230  (blas-level1-wrapx (iaxpy n alpha x incx offx y incy offy)
2231                     (y)
2232                     (lambda (i) (cond (conc "error code " i))))
2233 
2234  (blas-level1-sd-wrapx (dot n x incx y incy)
2235                        #f
2236                        (lambda (i) (cond (conc "error code " i))))
2237 
2238  (blas-level1-cz-wrapx (dotu n x incx y incy dotu)
2239                        #f
2240                        (lambda (i) (cond (conc "error code " i))))
2241 
2242  (blas-level1-cz-wrapx (dotc n x incx y incy dotc)
2243                        #f
2244                        (lambda (i) (cond (conc "error code " i))))
2245 
2246  (blas-level1-wrapx (nrm2 n x incx)
2247                     #f
2248                     (lambda (i) (cond (conc "error code " i))))
2249 
2250  (blas-level1-wrapx (asum n x incx)
2251                     #f
2252                     (lambda (i) (cond (conc "error code " i))))
2253 
2254  (blas-level1-wrapx (amax n x incx)
2255                     #f
2256                     (lambda (i) (cond (conc "error code " i))))
2257
2258
2259)
Note: See TracBrowser for help on using the repository browser.