H2Lib  3.0
blas.h
Go to the documentation of this file.
1 /* ------------------------------------------------------------
2  This is the file "blas.h" of the H2Lib package.
3  All rights reserved, Sven Christophersen 2015
4  ------------------------------------------------------------ */
5 
12 #ifndef BLAS_H_
13 #define BLAS_H_
14 
15 #ifdef USE_BLAS
16 
17 /* C STD LIBRARY */
18 /* CORE 0 */
19 #include "settings.h"
20 /* CORE 1 */
21 /* CORE 2 */
22 /* CORE 3 */
23 /* SIMPLE */
24 /* PARTICLES */
25 /* BEM */
26 
35 /****************************************************
36  * BLAS level 1
37  ****************************************************/
38 
39 /****************************************************
40  * _DOT / _RDOT
41  ****************************************************/
42 
55 #ifdef USE_FLOAT
56 #ifdef USE_COMPLEX
57 
58 IMPORT_PREFIX float _Complex
59 cdotc_(const unsigned *n, const float _Complex * x, const unsigned *incx,
60  const float _Complex * y, const unsigned *incy);
62 #define h2_dot(n, x, incx, y, incy) cdotc_(n, x, incx, y, incy)
63 #else
64 
65 IMPORT_PREFIX float
66 sdot_(const unsigned *n,
67  const float *x, const unsigned *incx, const float *y,
68  const unsigned *incy);
70 #define h2_dot(n, x, incx, y, incy) sdot_(n, x, incx, y, incy)
71 #endif
72 #else
73 #ifdef USE_COMPLEX
74 
75 IMPORT_PREFIX double _Complex
76 zdotc_(const unsigned *n, const double _Complex * x, const unsigned *incx,
77  const double _Complex * y, const unsigned *incy);
79 #define h2_dot(n, x, incx, y, incy) zdotc_(n, x, incx, y, incy)
80 #else
81 
82 IMPORT_PREFIX double
83 ddot_(const unsigned *n, const double *x, const unsigned *incx, const double *y,
84  const unsigned *incy);
86 #define h2_dot(n, x, incx, y, incy) ddot_(n, x, incx, y, incy)
87 #endif
88 #endif
89 
102 #ifdef USE_FLOAT
103 #ifdef USE_COMPLEX
104 
105 IMPORT_PREFIX float
106 sdot_(const unsigned *n, const float *x, const unsigned *incx, const float *y,
107  const unsigned *incy);
109 #define h2_rdot(n, x, incx, y, incy) sdot_(n, x, incx, y, incy)
110 #else
111 #define h2_rdot(n, x, incx, y, incy) sdot_(n, x, incx, y, incy)
112 #endif
113 #else
114 #ifdef USE_COMPLEX
115 
116 IMPORT_PREFIX double
117 ddot_(const unsigned *n, const double *x, const unsigned *incx, const double *y,
118  const unsigned *incy);
120 #define h2_rdot(n, x, incx, y, incy) ddot_(n, x, incx, y, incy)
121 #else
122 #define h2_rdot(n, x, incx, y, incy) ddot_(n, x, incx, y, incy)
123 #endif
124 #endif
125 
126 /****************************************************
127  * _AXPY /_RAXPY
128  ****************************************************/
129 
141 #ifdef USE_FLOAT
142 #ifdef USE_COMPLEX
143 
144 IMPORT_PREFIX void
145 caxpy_(const unsigned *n, const float _Complex * alpha,
146  const float _Complex * x, const unsigned *incx, float _Complex * y,
147  const unsigned *incy);
149 #define h2_axpy(n, alpha, x, incx, y, incy) caxpy_(n, alpha, x, incx, y, incy)
150 #else
151 
152 IMPORT_PREFIX void
153 saxpy_(const unsigned *n, const float * alpha,
154  const float * x, const unsigned *incx, float * y, const unsigned *incy);
156 #define h2_axpy(n, alpha, x, incx, y, incy) saxpy_(n, alpha, x, incx, y, incy)
157 #endif
158 #else
159 #ifdef USE_COMPLEX
160 
161 IMPORT_PREFIX void
162 zaxpy_(const unsigned *n, const double _Complex *alpha,
163  const double _Complex *x, const unsigned *incx, double _Complex *y,
164  const unsigned *incy);
166 #define h2_axpy(n, alpha, x, incx, y, incy) zaxpy_(n, alpha, x, incx, y, incy)
167 #else
168 
169 IMPORT_PREFIX void
170 daxpy_(const unsigned *n, const double *alpha, const double *x,
171  const unsigned *incx, double *y, const unsigned *incy);
173 #define h2_axpy(n, alpha, x, incx, y, incy) daxpy_(n, alpha, x, incx, y, incy)
174 #endif
175 #endif
176 
188 #ifdef USE_FLOAT
189 #ifdef USE_COMPLEX
190 
191 IMPORT_PREFIX void
192 saxpy_(const unsigned *n, const float * alpha, const float * x,
193  const unsigned *incx, float * y, const unsigned *incy);
195 #define h2_raxpy(n, alpha, x, incx, y, incy) saxpy_(n, alpha, x, incx, y, incy)
196 #else
197 #define h2_raxpy(n, alpha, x, incx, y, incy) saxpy_(n, alpha, x, incx, y, incy)
198 #endif
199 #else
200 #ifdef USE_COMPLEX
201 
202 IMPORT_PREFIX void
203 daxpy_(const unsigned *n, const double *alpha, const double *x,
204  const unsigned *incx, double *y, const unsigned *incy);
206 #define h2_raxpy(n, alpha, x, incx, y, incy) daxpy_(n, alpha, x, incx, y, incy)
207 #else
208 #define h2_raxpy(n, alpha, x, incx, y, incy) daxpy_(n, alpha, x, incx, y, incy)
209 #endif
210 #endif
211 
212 /****************************************************
213  * _SCAL / _RSCAL
214  ****************************************************/
215 
224 #ifdef USE_FLOAT
225 #ifdef USE_COMPLEX
226 
227 IMPORT_PREFIX void
228 cscal_(const unsigned *n, const float _Complex *alpha, float _Complex *x,
229  const unsigned *incx);
231 #define h2_scal(n, alpha, x, incx) cscal_(n, alpha, x, incx)
232 #else
233 
234 IMPORT_PREFIX void
235 sscal_(const unsigned *n, const float *alpha, float *x, const unsigned *incx);
237 #define h2_scal(n, alpha, x, incx) sscal_(n, alpha, x, incx)
238 #endif
239 #else
240 #ifdef USE_COMPLEX
241 
242 IMPORT_PREFIX void
243 zscal_(const unsigned *n, const double _Complex *alpha, double _Complex *x,
244  const unsigned *incx);
246 #define h2_scal(n, alpha, x, incx) zscal_(n, alpha, x, incx)
247 #else
248 
249 IMPORT_PREFIX void
250 dscal_(const unsigned *n, const double *alpha, double *x, const unsigned *incx);
252 #define h2_scal(n, alpha, x, incx) dscal_(n, alpha, x, incx)
253 #endif
254 #endif
255 
264 #ifdef USE_FLOAT
265 #ifdef USE_COMPLEX
266 
267 IMPORT_PREFIX void
268 csscal_(const unsigned *n, const float *alpha, float _Complex *x,
269  const unsigned *incx);
271 #define h2_rscal(n, alpha, x, incx) csscal_(n, alpha, x, incx)
272 #else
273 #define h2_rscal(n, alpha, x, incx) sscal_(n, alpha, x, incx)
274 #endif
275 #else
276 #ifdef USE_COMPLEX
277 
278 IMPORT_PREFIX void
279 zdscal_(const unsigned *n, const double *alpha, double _Complex *x,
280  const unsigned *incx);
282 #define h2_rscal(n, alpha, x, incx) zdscal_(n, alpha, x, incx)
283 #else
284 #define h2_rscal(n, alpha, x, incx) dscal_(n, alpha, x, incx)
285 #endif
286 #endif
287 
288 /****************************************************
289  * _NRM2
290  ****************************************************/
291 
300 #ifdef USE_FLOAT
301 #ifdef USE_COMPLEX
302 
303 IMPORT_PREFIX float
304 scnrm2_(const unsigned *n, const float _Complex *x, const unsigned *incx);
306 #define h2_nrm2(n, x, incx) scnrm2_(n, x, incx)
307 #else
308 
309 IMPORT_PREFIX float
310 snrm2_(const unsigned *n, const float *x, const unsigned *incx);
312 #define h2_nrm2(n, x, incx) snrm2_(n, x, incx)
313 #endif
314 #else
315 #ifdef USE_COMPLEX
316 
317 IMPORT_PREFIX double
318 dznrm2_(const unsigned *n, const double _Complex *x, const unsigned *incx);
320 #define h2_nrm2(n, x, incx) dznrm2_(n, x, incx)
321 #else
322 
323 IMPORT_PREFIX double
324 dnrm2_(const unsigned *n, const double *x, const unsigned *incx);
326 #define h2_nrm2(n, x, incx) dnrm2_(n, x, incx)
327 #endif
328 #endif
329 
330 /****************************************************
331  * _SWAP
332  ****************************************************/
333 
343 #ifdef USE_FLOAT
344 #ifdef USE_COMPLEX
345 
346 IMPORT_PREFIX void
347 cswap_(const unsigned *n, float _Complex *x, const unsigned *incx,
348  float _Complex *y, const unsigned *incy);
350 #define h2_swap(n, x, incx, y, incy) cswap_(n, x, incx, y, incy)
351 #else
352 
353 IMPORT_PREFIX void
354 sswap_(const unsigned *n, float *x, const unsigned *incx, float *y,
355  const unsigned *incy);
357 #define h2_swap(n, x, incx, y, incy) sswap_(n, x, incx, y, incy)
358 #endif
359 #else
360 #ifdef USE_COMPLEX
361 
362 IMPORT_PREFIX void
363 zswap_(const unsigned *n, double _Complex *x, const unsigned *incx,
364  double _Complex *y, const unsigned *incy);
366 #define h2_swap(n, x, incx, y, incy) zswap_(n, x, incx, y, incy)
367 #else
368 
369 IMPORT_PREFIX void
370 dswap_(const unsigned *n, double *x, const unsigned *incx, double *y,
371  const unsigned *incy);
373 #define h2_swap(n, x, incx, y, incy) dswap_(n, x, incx, y, incy)
374 #endif
375 #endif
376 
377 /****************************************************
378  * BLAS level 2
379  ****************************************************/
380 
381 /****************************************************
382  * _GEMV
383  ****************************************************/
384 
411 #ifdef USE_FLOAT
412 #ifdef USE_COMPLEX
413 
414 IMPORT_PREFIX void
415 cgemv_(const char *trans, const unsigned *m, const unsigned *n,
416  const float _Complex *alpha, const float _Complex *a, const unsigned *lda,
417  const float _Complex *x, const unsigned *incx, const float _Complex *beta,
418  float _Complex *y, const unsigned *incy);
420 #define h2_gemv(trans, m, n, alpha, a, lda, x, incx, beta, y, incy) cgemv_(trans, m, n, alpha, a, lda, x, incx, beta, y, incy)
421 #else
422 
423 IMPORT_PREFIX void
424 sgemv_(const char *trans, const unsigned *m, const unsigned *n,
425  const float *alpha, const float *a, const unsigned *lda, const float *x,
426  const unsigned *incx, const float *beta, float *y, const unsigned *incy);
428 #define h2_gemv(trans, m, n, alpha, a, lda, x, incx, beta, y, incy) sgemv_(trans, m, n, alpha, a, lda, x, incx, beta, y, incy)
429 #endif
430 #else
431 #ifdef USE_COMPLEX
432 
433 IMPORT_PREFIX void
434 zgemv_(const char *trans, const unsigned *m, const unsigned *n,
435  const double _Complex *alpha, const double _Complex *a, const unsigned *lda,
436  const double _Complex *x, const unsigned *incx, const double _Complex *beta,
437  double _Complex *y, const unsigned *incy);
439 #define h2_gemv(trans, m, n, alpha, a, lda, x, incx, beta, y, incy) zgemv_(trans, m, n, alpha, a, lda, x, incx, beta, y, incy)
440 #else
441 
442 IMPORT_PREFIX void
443 dgemv_(const char *trans, const unsigned *m, const unsigned *n,
444  const double *alpha, const double *a, const unsigned *lda, const double *x,
445  const unsigned *incx, const double *beta, double *y, const unsigned *incy);
447 #define h2_gemv(trans, m, n, alpha, a, lda, x, incx, beta, y, incy) dgemv_(trans, m, n, alpha, a, lda, x, incx, beta, y, incy)
448 #endif
449 #endif
450 
451 /****************************************************
452  * _TRMV
453  ****************************************************/
454 
479 #ifdef USE_FLOAT
480 #ifdef USE_COMPLEX
481 
482 IMPORT_PREFIX void
483 ctrmv_(const char *uplo, const char *trans, const char *diag, const unsigned *n,
484  const float _Complex *a, const unsigned *lda, float _Complex *x,
485  const unsigned *incx);
487 #define h2_trmv(uplo, trans, diag, n, a, lda, x, incx) ctrmv_(uplo, trans, diag, n, a, lda, x, incx)
488 #else
489 
490 IMPORT_PREFIX void
491 strmv_(const char *uplo, const char *trans, const char *diag, const unsigned *n,
492  const float *a, const unsigned *lda, float *x, const unsigned *incx);
494 #define h2_trmv(uplo, trans, diag, n, a, lda, x, incx) strmv_(uplo, trans, diag, n, a, lda, x, incx)
495 #endif
496 #else
497 #ifdef USE_COMPLEX
498 
499 IMPORT_PREFIX void
500 ztrmv_(const char *uplo, const char *trans, const char *diag, const unsigned *n,
501  const double _Complex *a, const unsigned *lda, double _Complex *x,
502  const unsigned *incx);
504 #define h2_trmv(uplo, trans, diag, n, a, lda, x, incx) ztrmv_(uplo, trans, diag, n, a, lda, x, incx)
505 #else
506 
507 IMPORT_PREFIX void
508 dtrmv_(const char *uplo, const char *trans, const char *diag, const unsigned *n,
509  const double *a, const unsigned *lda, double *x, const unsigned *incx);
511 #define h2_trmv(uplo, trans, diag, n, a, lda, x, incx) dtrmv_(uplo, trans, diag, n, a, lda, x, incx)
512 #endif
513 #endif
514 
515 /****************************************************
516  * _GER
517  ****************************************************/
518 
534 #ifdef USE_FLOAT
535 #ifdef USE_COMPLEX
536 
537 IMPORT_PREFIX void
538 cgerc_(const unsigned *m, const unsigned *n, const float _Complex *alpha,
539  const float _Complex *x, const unsigned *incx, const float _Complex *y,
540  const unsigned *incy, float _Complex *a, const unsigned *lda);
542 #define h2_ger(m, n, alpha, x, incx, y, incy, a, lda) cgerc_(m, n, alpha, x, incx, y, incy, a, lda)
543 #else
544 
545 IMPORT_PREFIX void
546 sger_(const unsigned *m, const unsigned *n, const float *alpha,
547  const float *x, const unsigned *incx, const float *y,
548  const unsigned *incy, float *a, const unsigned *lda);
550 #define h2_ger(m, n, alpha, x, incx, y, incy, a, lda) sger_(m, n, alpha, x, incx, y, incy, a, lda)
551 #endif
552 #else
553 #ifdef USE_COMPLEX
554 
555 IMPORT_PREFIX void
556 zgerc_(const unsigned *m, const unsigned *n, const double _Complex *alpha,
557  const double _Complex *x, const unsigned *incx, const double _Complex *y,
558  const unsigned *incy, double _Complex *a, const unsigned *lda);
560 #define h2_ger(m, n, alpha, x, incx, y, incy, a, lda) zgerc_(m, n, alpha, x, incx, y, incy, a, lda)
561 #else
562 
563 IMPORT_PREFIX void
564 dger_(const unsigned *m, const unsigned *n, const double *alpha,
565  const double *x, const unsigned *incx, const double *y,
566  const unsigned *incy, double *a, const unsigned *lda);
568 #define h2_ger(m, n, alpha, x, incx, y, incy, a, lda) dger_(m, n, alpha, x, incx, y, incy, a, lda)
569 #endif
570 #endif
571 
587 #ifdef USE_FLOAT
588 #ifdef USE_COMPLEX
589 #define h2_gerc(m, n, alpha, x, incx, y, incy, a, lda) cgerc_(m, n, alpha, x, incx, y, incy, a, lda)
590 #else
591 #define h2_gerc(m, n, alpha, x, incx, y, incy, a, lda) sger_(m, n, alpha, x, incx, y, incy, a, lda)
592 #endif
593 #else
594 #ifdef USE_COMPLEX
595 #define h2_gerc(m, n, alpha, x, incx, y, incy, a, lda) zgerc_(m, n, alpha, x, incx, y, incy, a, lda)
596 #else
597 #define h2_gerc(m, n, alpha, x, incx, y, incy, a, lda) dger_(m, n, alpha, x, incx, y, incy, a, lda)
598 #endif
599 #endif
600 
616 #ifdef USE_FLOAT
617 #ifdef USE_COMPLEX
618 
619 IMPORT_PREFIX void
620 cgeru_(const unsigned *m, const unsigned *n, const float _Complex *alpha,
621  const float _Complex *x, const unsigned *incx, const float _Complex *y,
622  const unsigned *incy, float _Complex *a, const unsigned *lda);
624 #define h2_geru(m, n, alpha, x, incx, y, incy, a, lda) cgeru_(m, n, alpha, x, incx, y, incy, a, lda)
625 #else
626 #define h2_geru(m, n, alpha, x, incx, y, incy, a, lda) sger_(m, n, alpha, x, incx, y, incy, a, lda)
627 #endif
628 #else
629 #ifdef USE_COMPLEX
630 
631 IMPORT_PREFIX void
632 zgeru_(const unsigned *m, const unsigned *n, const double _Complex *alpha,
633  const double _Complex *x, const unsigned *incx, const double _Complex *y,
634  const unsigned *incy, double _Complex *a, const unsigned *lda);
636 #define h2_geru(m, n, alpha, x, incx, y, incy, a, lda) zgeru_(m, n, alpha, x, incx, y, incy, a, lda)
637 #else
638 #define h2_geru(m, n, alpha, x, incx, y, incy, a, lda) dger_(m, n, alpha, x, incx, y, incy, a, lda)
639 #endif
640 #endif
641 
642 /****************************************************
643  * _SYR
644  ****************************************************/
645 
661 #ifdef USE_FLOAT
662 #ifdef USE_COMPLEX
663 
664 IMPORT_PREFIX void
665 cher_(const char *uplo, const unsigned *n, const float *alpha,
666  const float _Complex *x, const unsigned *incx, float _Complex *a,
667  const unsigned *lda);
669 #define h2_syr(uplo, n, alpha, x, incx, a, lda) cher_(uplo, n, alpha, x, incx, a, lda)
670 #else
671 
672 IMPORT_PREFIX void
673 ssyr_(const char *uplo, const unsigned *n, const float *alpha, const float *x,
674  const unsigned *incx, float *a, const unsigned *lda);
676 #define h2_syr(uplo, n, alpha, x, incx, a, lda) ssyr_(uplo, n, alpha, x, incx, a, lda)
677 #endif
678 #else
679 #ifdef USE_COMPLEX
680 
681 IMPORT_PREFIX void
682 zher_(const char *uplo, const unsigned *n, const double *alpha,
683  const double _Complex *x, const unsigned *incx, double _Complex *a,
684  const unsigned *lda);
686 #define h2_syr(uplo, n, alpha, x, incx, a, lda) zher_(uplo, n, alpha, x, incx, a, lda)
687 #else
688 
689 IMPORT_PREFIX void
690 dsyr_(const char *uplo, const unsigned *n, const double *alpha, const double *x,
691  const unsigned *incx, double *a, const unsigned *lda);
693 #define h2_syr(uplo, n, alpha, x, incx, a, lda) dsyr_(uplo, n, alpha, x, incx, a, lda)
694 #endif
695 #endif
696 
697 /****************************************************
698  * BLAS level 3
699  ****************************************************/
700 
701 /****************************************************
702  * _GEMM
703  ****************************************************/
704 
734 #ifdef USE_FLOAT
735 #ifdef USE_COMPLEX
736 
737 IMPORT_PREFIX void
738 cgemm_(const char *transa, const char *transb, const unsigned *m,
739  const unsigned *n, const unsigned *k, const float _Complex *alpha,
740  const float _Complex *a, const unsigned *lda, const float _Complex *b,
741  const unsigned *ldb, const float _Complex *beta, float _Complex *c,
742  const unsigned *ldc);
744 #define h2_gemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc) cgemm_(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
745 #else
746 
747 IMPORT_PREFIX void
748 sgemm_(const char *transa, const char *transb, const unsigned *m,
749  const unsigned *n, const unsigned *k, const float *alpha, const float *a,
750  const unsigned *lda, const float *b, const unsigned *ldb,
751  const float *beta, float *c, const unsigned *ldc);
753 #define h2_gemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc) sgemm_(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
754 #endif
755 #else
756 #ifdef USE_COMPLEX
757 
758 IMPORT_PREFIX void
759 zgemm_(const char *transa, const char *transb, const unsigned *m,
760  const unsigned *n, const unsigned *k, const double _Complex *alpha,
761  const double _Complex *a, const unsigned *lda, const double _Complex *b,
762  const unsigned *ldb, const double _Complex *beta, double _Complex *c,
763  const unsigned *ldc);
765 #define h2_gemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc) zgemm_(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
766 #else
767 
768 IMPORT_PREFIX void
769 dgemm_(const char *transa, const char *transb, const unsigned *m,
770  const unsigned *n, const unsigned *k, const double *alpha, const double *a,
771  const unsigned *lda, const double *b, const unsigned *ldb,
772  const double *beta, double *c, const unsigned *ldc);
774 #define h2_gemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc) dgemm_(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
775 #endif
776 #endif
777 
778 /****************************************************
779  * _TRMM
780  ****************************************************/
781 
813 #ifdef USE_FLOAT
814 #ifdef USE_COMPLEX
815 
816 IMPORT_PREFIX void
817 ctrmm_(const char *side, const char *uplo, const char *transa, const char *diag,
818  const unsigned *m, const unsigned *n, const float _Complex *alpha,
819  const float _Complex *a, const unsigned *lda, float _Complex *b,
820  const unsigned *ldb);
822 #define h2_trmm(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb) ctrmm_(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb)
823 #else
824 
825 IMPORT_PREFIX void
826 strmm_(const char *side, const char *uplo, const char *transa, const char *diag,
827  const unsigned *m, const unsigned *n, const float *alpha, const float *a,
828  const unsigned *lda, float *b, const unsigned *ldb);
830 #define h2_trmm(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb) strmm_(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb)
831 #endif
832 #else
833 #ifdef USE_COMPLEX
834 
835 IMPORT_PREFIX void
836 ztrmm_(const char *side, const char *uplo, const char *transa, const char *diag,
837  const unsigned *m, const unsigned *n, const double _Complex *alpha,
838  const double _Complex *a, const unsigned *lda, double _Complex *b,
839  const unsigned *ldb);
841 #define h2_trmm(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb) ztrmm_(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb)
842 #else
843 
844 IMPORT_PREFIX void
845 dtrmm_(const char *side, const char *uplo, const char *transa, const char *diag,
846  const unsigned *m, const unsigned *n, const double *alpha, const double *a,
847  const unsigned *lda, double *b, const unsigned *ldb);
849 #define h2_trmm(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb) dtrmm_(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb)
850 #endif
851 #endif
852 
853 /****************************************************
854  * _TRSM
855  ****************************************************/
856 
885 #ifdef USE_FLOAT
886 #ifdef USE_COMPLEX
887 
888 IMPORT_PREFIX void
889 ctrsm_(const char *side, const char *uplo, const char *transa, const char *diag,
890  const unsigned *m, const unsigned *n, const float _Complex *alpha,
891  const float _Complex *a, const unsigned *lda, _Complex float *b,
892  const unsigned *ldb);
894 #define h2_trsm(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb) ctrsm_(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb)
895 #else
896 
897 IMPORT_PREFIX void
898 strsm_(const char *side, const char *uplo, const char *transa, const char *diag,
899  const unsigned *m, const unsigned *n, const float *alpha, const float *a,
900  const unsigned *lda, float *b, const unsigned *ldb);
902 #define h2_trsm(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb) strsm_(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb)
903 #endif
904 #else
905 #ifdef USE_COMPLEX
906 
907 IMPORT_PREFIX void
908 ztrsm_(const char *side, const char *uplo, const char *transa, const char *diag,
909  const unsigned *m, const unsigned *n, const double _Complex *alpha,
910  const double _Complex *a, const unsigned *lda, _Complex double *b,
911  const unsigned *ldb);
913 #define h2_trsm(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb) ztrsm_(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb)
914 #else
915 
916 IMPORT_PREFIX void
917 dtrsm_(const char *side, const char *uplo, const char *transa, const char *diag,
918  const unsigned *m, const unsigned *n, const double *alpha, const double *a,
919  const unsigned *lda, double *b, const unsigned *ldb);
921 #define h2_trsm(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb) dtrsm_(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb)
922 #endif
923 #endif
924 
925 /****************************************************
926  * LAPACK
927  ****************************************************/
928 
929 /****************************************************
930  * _LACGV
931  ****************************************************/
932 
942 #ifdef USE_FLOAT
943 #ifdef USE_COMPLEX
944 
945 IMPORT_PREFIX void
946 clacgv_(const unsigned *n, float _Complex *x, const unsigned *incx);
948 #define h2_lacgv(n, x, incx) clacgv_(n, x, incx)
949 #else
950 #define h2_lacgv(n, x, incx) (void)
951 #endif
952 #else
953 #ifdef USE_COMPLEX
954 
955 IMPORT_PREFIX void
956 zlacgv_(const unsigned *n, double _Complex *x, const unsigned *incx);
958 #define h2_lacgv(n, x, incx) zlacgv_(n, x, incx)
959 #else
960 #define h2_lacgv(n, x, incx) (void)
961 #endif
962 #endif
963 
964 /****************************************************
965  * _POTRF
966  ****************************************************/
967 
987 #ifdef USE_FLOAT
988 #ifdef USE_COMPLEX
989 
990 IMPORT_PREFIX void
991 cpotrf_(const char *uplo, const unsigned *n, float _Complex *a,
992  const unsigned *lda, int *info);
994 #define h2_potrf(uplo, n, a, lda, info) cpotrf_(uplo, n, a, lda, info)
995 #else
996 
997 IMPORT_PREFIX void
998 spotrf_(const char *uplo, const unsigned *n, float *a, const unsigned *lda,
999  int *info);
1001 #define h2_potrf(uplo, n, a, lda, info) spotrf_(uplo, n, a, lda, info)
1002 #endif
1003 #else
1004 #ifdef USE_COMPLEX
1005 
1006 IMPORT_PREFIX void
1007 zpotrf_(const char *uplo, const unsigned *n, double _Complex *a,
1008  const unsigned *lda, int *info);
1010 #define h2_potrf(uplo, n, a, lda, info) zpotrf_(uplo, n, a, lda, info)
1011 #else
1012 
1013 IMPORT_PREFIX void
1014 dpotrf_(const char *uplo, const unsigned *n, double *a, const unsigned *lda,
1015  int *info);
1017 #define h2_potrf(uplo, n, a, lda, info) dpotrf_(uplo, n, a, lda, info)
1018 #endif
1019 #endif
1020 
1021 /****************************************************
1022  * _LARF / _LARFG
1023  ****************************************************/
1024 
1045 #ifdef USE_FLOAT
1046 #ifdef USE_COMPLEX
1047 
1048 IMPORT_PREFIX void
1049 clarf_(const char *side, const unsigned *m, const unsigned *n,
1050  const float _Complex *v, const unsigned *incv, const float _Complex *tau,
1051  float _Complex *c, const unsigned *ldc, float _Complex *work);
1053 #define h2_larf(side, m, n, v, incv, tau, c, ldc, work) clarf_(side, m, n, v, incv, tau, c, ldc, work)
1054 #else
1055 
1056 IMPORT_PREFIX void
1057 slarf_(const char *side, const unsigned *m, const unsigned *n, const float *v,
1058  const unsigned *incv, const float *tau, float *c, const unsigned *ldc,
1059  float *work);
1061 #define h2_larf(side, m, n, v, incv, tau, c, ldc, work) slarf_(side, m, n, v, incv, tau, c, ldc, work)
1062 #endif
1063 #else
1064 #ifdef USE_COMPLEX
1065 
1066 IMPORT_PREFIX void
1067 zlarf_(const char *side, const unsigned *m, const unsigned *n,
1068  const double _Complex *v, const unsigned *incv, const double _Complex *tau,
1069  double _Complex *c, const unsigned *ldc, double _Complex *work);
1071 #define h2_larf(side, m, n, v, incv, tau, c, ldc, work) zlarf_(side, m, n, v, incv, tau, c, ldc, work)
1072 #else
1073 
1074 IMPORT_PREFIX void
1075 dlarf_(const char *side, const unsigned *m, const unsigned *n, const double *v,
1076  const unsigned *incv, const double *tau, double *c, const unsigned *ldc,
1077  double *work);
1079 #define h2_larf(side, m, n, v, incv, tau, c, ldc, work) dlarf_(side, m, n, v, incv, tau, c, ldc, work)
1080 #endif
1081 #endif
1082 
1103 #ifdef USE_FLOAT
1104 #ifdef USE_COMPLEX
1105 
1106 IMPORT_PREFIX void
1107 clarfg_(const unsigned *n, float _Complex *alpha, float _Complex *x,
1108  const unsigned *incx, float _Complex *tau);
1110 #define h2_larfg(n, alpha, x, incx, tau) clarfg_(n, alpha, x, incx, tau)
1111 #else
1112 
1113 IMPORT_PREFIX void
1114 slarfg_(const unsigned *n, float *alpha, float *x, const unsigned *incx, float *tau);
1116 #define h2_larfg(n, alpha, x, incx, tau) slarfg_(n, alpha, x, incx, tau)
1117 #endif
1118 #else
1119 #ifdef USE_COMPLEX
1120 
1121 IMPORT_PREFIX void
1122 zlarfg_(const unsigned *n, double _Complex *alpha, double _Complex *x,
1123  const unsigned *incx, double _Complex *tau);
1125 #define h2_larfg(n, alpha, x, incx, tau) zlarfg_(n, alpha, x, incx, tau)
1126 #else
1127 
1128 IMPORT_PREFIX void
1129 dlarfg_(const unsigned *n, double *alpha, double *x, const unsigned *incx,
1130  double *tau);
1132 #define h2_larfg(n, alpha, x, incx, tau) dlarfg_(n, alpha, x, incx, tau)
1133 #endif
1134 #endif
1135 
1136 /****************************************************
1137  * _GEQRF
1138  ****************************************************/
1139 
1160 #ifdef USE_FLOAT
1161 #ifdef USE_COMPLEX
1162 
1163 IMPORT_PREFIX void
1164 cgeqrf_(const unsigned *m, const unsigned *n, float _Complex *a,
1165  const unsigned *lda, float _Complex *tau, float _Complex *work,
1166  const int *lwork, int *info);
1168 #define h2_geqrf(m, n, a, lda, tau, work, lwork, info) cgeqrf_(m, n, a, lda, tau, work, lwork, info)
1169 #else
1170 
1171 IMPORT_PREFIX void
1172 sgeqrf_(const unsigned *m, const unsigned *n, float *a, const unsigned *lda,
1173  float *tau, float *work, const int *lwork, int *info);
1175 #define h2_geqrf(m, n, a, lda, tau, work, lwork, info) sgeqrf_(m, n, a, lda, tau, work, lwork, info)
1176 #endif
1177 #else
1178 #ifdef USE_COMPLEX
1179 
1180 IMPORT_PREFIX void
1181 zgeqrf_(const unsigned *m, const unsigned *n, double _Complex *a,
1182  const unsigned *lda, double _Complex *tau, double _Complex *work,
1183  const int *lwork, int *info);
1185 #define h2_geqrf(m, n, a, lda, tau, work, lwork, info) zgeqrf_(m, n, a, lda, tau, work, lwork, info)
1186 #else
1187 
1188 IMPORT_PREFIX void
1189 dgeqrf_(const unsigned *m, const unsigned *n, double *a, const unsigned *lda,
1190  double *tau, double *work, const int *lwork, int *info);
1192 #define h2_geqrf(m, n, a, lda, tau, work, lwork, info) dgeqrf_(m, n, a, lda, tau, work, lwork, info)
1193 #endif
1194 #endif
1195 
1196 /****************************************************
1197  * _ORMQR
1198  ****************************************************/
1199 
1233 #ifdef USE_FLOAT
1234 #ifdef USE_COMPLEX
1235 
1236 IMPORT_PREFIX void
1237 cunmqr_(const char *side, const char *trans, const unsigned *m,
1238  const unsigned *n, const unsigned *k, const float _Complex *a,
1239  const unsigned *lda, const float _Complex *tau, float _Complex *c,
1240  const unsigned *ldc, float _Complex *work, const int *lwork, int *info);
1242 #define h2_ormqr(side, trans, m, n, k, a, lda, tau, c, ldc, work, lwork, info) cunmqr_(side, trans, m, n, k, a, lda, tau, c, ldc, work, lwork, info)
1243 #else
1244 
1245 IMPORT_PREFIX void
1246 sormqr_(const char *side, const char *trans, const unsigned *m,
1247  const unsigned *n, const unsigned *k, const float *a, const unsigned *lda,
1248  const float *tau, float *c, const unsigned *ldc, float *work,
1249  const int *lwork, int *info);
1251 #define h2_ormqr(side, trans, m, n, k, a, lda, tau, c, ldc, work, lwork, info) sormqr_(side, trans, m, n, k, a, lda, tau, c, ldc, work, lwork, info)
1252 #endif
1253 #else
1254 #ifdef USE_COMPLEX
1255 
1256 IMPORT_PREFIX void
1257 zunmqr_(const char *side, const char *trans, const unsigned *m,
1258  const unsigned *n, const unsigned *k, const double _Complex *a,
1259  const unsigned *lda, const double _Complex *tau, double _Complex *c,
1260  const unsigned *ldc, double _Complex *work, const int *lwork, int *info);
1262 #define h2_ormqr(side, trans, m, n, k, a, lda, tau, c, ldc, work, lwork, info) zunmqr_(side, trans, m, n, k, a, lda, tau, c, ldc, work, lwork, info)
1263 #else
1264 
1265 IMPORT_PREFIX void
1266 dormqr_(const char *side, const char *trans, const unsigned *m,
1267  const unsigned *n, const unsigned *k, const double *a, const unsigned *lda,
1268  const double *tau, double *c, const unsigned *ldc, double *work,
1269  const int *lwork, int *info);
1271 #define h2_ormqr(side, trans, m, n, k, a, lda, tau, c, ldc, work, lwork, info) dormqr_(side, trans, m, n, k, a, lda, tau, c, ldc, work, lwork, info)
1272 #endif
1273 #endif
1274 
1275 /****************************************************
1276  * _ORGQR
1277  ****************************************************/
1278 
1299 #ifdef USE_FLOAT
1300 #ifdef USE_COMPLEX
1301 
1302 IMPORT_PREFIX void
1303 cungqr_(const unsigned *m, const unsigned *n, const unsigned *k,
1304  float _Complex *a, const unsigned *lda, float _Complex *tau,
1305  float _Complex *work, const unsigned *lwork, unsigned *info);
1307 #define h2_orgqr(m, n, k, a, lda, tau, work, lwork, info) cungqr_(m, n, k, a, lda, tau, work, lwork, info)
1308 #else
1309 
1310 IMPORT_PREFIX void
1311 sorgqr_(const unsigned *m, const unsigned *n, const unsigned *k, float *a,
1312  const unsigned *lda, float *tau, float *work, const unsigned *lwork,
1313  unsigned *info);
1315 #define h2_orgqr(m, n, k, a, lda, tau, work, lwork, info) sorgqr_(m, n, k, a, lda, tau, work, lwork, info)
1316 #endif
1317 #else
1318 #ifdef USE_COMPLEX
1319 
1320 IMPORT_PREFIX void
1321 zungqr_(const unsigned *m, const unsigned *n, const unsigned *k,
1322  double _Complex *a, const unsigned *lda, double _Complex *tau,
1323  double _Complex *work, const unsigned *lwork, unsigned *info);
1325 #define h2_orgqr(m, n, k, a, lda, tau, work, lwork, info) zungqr_(m, n, k, a, lda, tau, work, lwork, info)
1326 #else
1327 
1328 IMPORT_PREFIX void
1329 dorgqr_(const unsigned *m, const unsigned *n, const unsigned *k, double *a,
1330  const unsigned *lda, double *tau, double *work, const unsigned *lwork,
1331  unsigned *info);
1333 #define h2_orgqr(m, n, k, a, lda, tau, work, lwork, info) dorgqr_(m, n, k, a, lda, tau, work, lwork, info)
1334 #endif
1335 #endif
1336 
1337 /****************************************************
1338  * _STEQR
1339  ****************************************************/
1340 
1359 #ifdef USE_FLOAT
1360 #ifdef USE_COMPLEX
1361 
1362 IMPORT_PREFIX void
1363 csteqr_(const char *compz, const unsigned *n, float *d, float *e,
1364  float _Complex *z, const unsigned *ldz, float *work, int *info);
1366 #define h2_steqr(compz, n, d, e, z, ldz, work, info) csteqr_(compz, n, d, e, z, ldz, work, info)
1367 #else
1368 
1369 IMPORT_PREFIX void
1370 ssteqr_(const char *compz,const unsigned *n,float *d,float *e, float *z,
1371  const unsigned *ldz, float *work, int *info);
1373 #define h2_steqr(compz, n, d, e, z, ldz, work, info) ssteqr_(compz, n, d, e, z, ldz, work, info)
1374 #endif
1375 #else
1376 #ifdef USE_COMPLEX
1377 
1378 IMPORT_PREFIX void
1379 zsteqr_(const char *compz, const unsigned *n, double *d, double *e,
1380  double _Complex *z, const unsigned *ldz, double *work, int *info);
1382 #define h2_steqr(compz, n, d, e, z, ldz, work, info) zsteqr_(compz, n, d, e, z, ldz, work, info)
1383 #else
1384 
1385 IMPORT_PREFIX void
1386 dsteqr_(const char *compz, const unsigned *n, double *d, double *e, double *z,
1387  const unsigned *ldz, double *work, int *info);
1389 #define h2_steqr(compz, n, d, e, z, ldz, work, info) dsteqr_(compz, n, d, e, z, ldz, work, info)
1390 #endif
1391 #endif
1392 
1393 /****************************************************
1394  * _STEV
1395  ****************************************************/
1396 
1413 #ifdef USE_FLOAT
1414 
1415 IMPORT_PREFIX void
1416 sstev_(const char *jobz, const unsigned *n, float *d, float *e, float *z,
1417  const unsigned *ldz, float *work, int *info);
1419 #define h2_stev(jobz, n, d, e, z, ldz, work, info) sstev_(jobz, n, d, e, z, ldz, work, info)
1420 #else
1421 
1422 IMPORT_PREFIX void
1423 dstev_(const char *jobz, const unsigned *n, double *d, double *e, double *z,
1424  const unsigned *ldz, double *work, int *info);
1426 #define h2_stev(jobz, n, d, e, z, ldz, work, info) dstev_(jobz, n, d, e, z, ldz, work, info)
1427 #endif
1428 
1429 /****************************************************
1430  * _GESVD
1431  ****************************************************/
1432 
1466 #ifdef USE_FLOAT
1467 #ifdef USE_COMPLEX
1468 
1469 IMPORT_PREFIX void
1470 cgesvd_(const char *jobu, const char *jobvt, const unsigned *m,
1471  const unsigned *n, float _Complex *a, const unsigned *lda, float *s,
1472  float _Complex *u, const unsigned *ldu, float _Complex *vt,
1473  const unsigned *ldvt, float _Complex *work, const unsigned *lwork,
1474  float *rwork, int *info);
1476 #define h2_gesvd(jobu, jobvt, m, n, a, lda, s, u, ldu, vt, ldvt, work, lwork, rwork, info) cgesvd_(jobu, jobvt, m, n, a, lda, s, u, ldu, vt, ldvt, work, lwork, rwork, info)
1477 #else
1478 
1479 IMPORT_PREFIX void
1480 sgesvd_(const char *jobu, const char *jobvt, const unsigned *m,
1481  const unsigned *n, float *a, const unsigned *lda, float *s, float *u,
1482  const unsigned *ldu, float *vt, const unsigned *ldvt, float *work,
1483  const unsigned *lwork, int *info);
1485 #define h2_gesvd(jobu, jobvt, m, n, a, lda, s, u, ldu, vt, ldvt, work, lwork, info) sgesvd_(jobu, jobvt, m, n, a, lda, s, u, ldu, vt, ldvt, work, lwork, info)
1486 #endif
1487 #else
1488 #ifdef USE_COMPLEX
1489 
1490 IMPORT_PREFIX void
1491 zgesvd_(const char *jobu, const char *jobvt, const unsigned *m,
1492  const unsigned *n, double _Complex *a, const unsigned *lda, double *s,
1493  double _Complex *u, const unsigned *ldu, double _Complex *vt,
1494  const unsigned *ldvt, double _Complex *work, const unsigned *lwork,
1495  double *rwork, int *info);
1497 #define h2_gesvd(jobu, jobvt, m, n, a, lda, s, u, ldu, vt, ldvt, work, lwork, rwork, info) zgesvd_(jobu, jobvt, m, n, a, lda, s, u, ldu, vt, ldvt, work, lwork, rwork, info)
1498 #else
1499 
1500 IMPORT_PREFIX void
1501 dgesvd_(const char *jobu, const char *jobvt, const unsigned *m,
1502  const unsigned *n, double *a, const unsigned *lda, double *s, double *u,
1503  const unsigned *ldu, double *vt, const unsigned *ldvt, double *work,
1504  const unsigned *lwork, int *info);
1506 #define h2_gesvd(jobu, jobvt, m, n, a, lda, s, u, ldu, vt, ldvt, work, lwork, info) dgesvd_(jobu, jobvt, m, n, a, lda, s, u, ldu, vt, ldvt, work, lwork, info)
1507 #endif
1508 #endif
1509 
1510 /****************************************************
1511  * _DBSQR
1512  ****************************************************/
1513 
1549 #ifdef USE_FLOAT
1550 #ifdef USE_COMPLEX
1551 
1552 IMPORT_PREFIX void
1553 cbdsqr_(const char *uplo, const unsigned *n, const unsigned *ncvt,
1554  const unsigned *nru, const unsigned *ncc, float *d, float *e,
1555  float _Complex *vt, const unsigned *ldvt, float _Complex *u,
1556  const unsigned *ldu, float _Complex *c, const unsigned *ldc, float *rwork,
1557  int *info);
1559 #define h2_bdsqr(uplo, n, ncvt, nru, ncc, d, e, vt, ldvt, u, ldu, c, ldc, rwork, info) cbdsqr_(uplo, n, ncvt, nru, ncc, d, e, vt, ldvt, u, ldu, c, ldc, rwork, info)
1560 #else
1561 
1562 IMPORT_PREFIX void
1563 sbdsqr_(const char *uplo, const unsigned *n, const unsigned *ncvt,
1564  const unsigned *nru, const unsigned *ncc, float *d, float *e, float *vt,
1565  const unsigned *ldvt, float *u, const unsigned *ldu, float *c,
1566  const unsigned *ldc, float *rwork, int *info);
1568 #define h2_bdsqr(uplo, n, ncvt, nru, ncc, d, e, vt, ldvt, u, ldu, c, ldc, rwork, info) sbdsqr_(uplo, n, ncvt, nru, ncc, d, e, vt, ldvt, u, ldu, c, ldc, rwork, info)
1569 #endif
1570 #else
1571 #ifdef USE_COMPLEX
1572 
1573 IMPORT_PREFIX void
1574 zbdsqr_(const char *uplo, const unsigned *n, const unsigned *ncvt,
1575  const unsigned *nru, const unsigned *ncc, double *d, double *e,
1576  double _Complex *vt, const unsigned *ldvt, double _Complex *u,
1577  const unsigned *ldu, double _Complex *c, const unsigned *ldc, double *rwork,
1578  int *info);
1580 #define h2_bdsqr(uplo, n, ncvt, nru, ncc, d, e, vt, ldvt, u, ldu, c, ldc, rwork, info) zbdsqr_(uplo, n, ncvt, nru, ncc, d, e, vt, ldvt, u, ldu, c, ldc, rwork, info)
1581 #else
1582 
1583 IMPORT_PREFIX void
1584 dbdsqr_(const char *uplo, const unsigned *n, const unsigned *ncvt,
1585  const unsigned *nru, const unsigned *ncc, double *d, double *e, double *vt,
1586  const unsigned *ldvt, double *u, const unsigned *ldu, double *c,
1587  const unsigned *ldc, double *rwork, int *info);
1589 #define h2_bdsqr(uplo, n, ncvt, nru, ncc, d, e, vt, ldvt, u, ldu, c, ldc, rwork, info) dbdsqr_(uplo, n, ncvt, nru, ncc, d, e, vt, ldvt, u, ldu, c, ldc, rwork, info)
1590 #endif
1591 #endif
1592 
1593 /****************************************************
1594  * _SYEV / _HEEV
1595  ****************************************************/
1596 
1619 #ifdef USE_FLOAT
1620 #ifdef USE_COMPLEX
1621 
1622 IMPORT_PREFIX void
1623 cheev_(const char *jobz, const char *uplo, const unsigned *n, float _Complex *a,
1624  const unsigned *lda, const float *w, float _Complex *work,
1625  const unsigned *lwork, float *rwork, int *info);
1627 #define h2_heev(jobz, uplo, n, a, lda, w, work, lwork, rwork, info) cheev_(jobz, uplo, n, a, lda, w, work, lwork, rwork, info)
1628 #else
1629 
1630 IMPORT_PREFIX void
1631 ssyev_(const char *jobz, const char *uplo, const unsigned *n, float *a,
1632  const unsigned *lda, const float *w, float *work, const unsigned *lwork,
1633  int *info);
1635 #define h2_heev(jobz, uplo, n, a, lda, w, work, lwork, rwork, info) ssyev_(jobz, uplo, n, a, lda, w, work, lwork, info)
1636 #endif
1637 #else
1638 #ifdef USE_COMPLEX
1639 
1640 IMPORT_PREFIX void
1641 zheev_(const char *jobz, const char *uplo, const unsigned *n,
1642  double _Complex *a, const unsigned *lda, const double *w,
1643  double _Complex *work, const unsigned *lwork, double *rwork, int *info);
1645 #define h2_heev(jobz, uplo, n, a, lda, w, work, lwork, rwork, info) zheev_(jobz, uplo, n, a, lda, w, work, lwork, rwork, info)
1646 #else
1647 
1648 IMPORT_PREFIX void
1649 dsyev_(const char *jobz, const char *uplo, const unsigned *n, double *a,
1650  const unsigned *lda, const double *w, double *work, const unsigned *lwork,
1651  int *info);
1653 #define h2_heev(jobz, uplo, n, a, lda, w, work, lwork, rwork, info) dsyev_(jobz, uplo, n, a, lda, w, work, lwork, info)
1654 #endif
1655 #endif
1656 
1657 #endif
1658 
1663 #endif
1664 
#define IMPORT_PREFIX
Prefix for external declarations.
Definition: settings.h:50