Teuchos - Trilinos Tools Package Version of the Day
Loading...
Searching...
No Matches
Teuchos_BLAS_wrappers.hpp
Go to the documentation of this file.
1// @HEADER
2// *****************************************************************************
3// Teuchos: Common Tools Package
4//
5// Copyright 2004 NTESS and the Teuchos contributors.
6// SPDX-License-Identifier: BSD-3-Clause
7// *****************************************************************************
8// @HEADER
9
10#ifndef _TEUCHOS_BLAS_WRAPPERS_HPP_
11#define _TEUCHOS_BLAS_WRAPPERS_HPP_
12
14#ifdef _MSC_VER
15/* disable warning for C-linkage returning complex class */
16#pragma warning ( disable : 4190 )
17#endif
18
24/* A) Define PREFIX and Teuchos_fcd based on platform. */
25
26#if defined(INTEL_CXML)
27# define PREFIX __stdcall
28# define Teuchos_fcd const char *, unsigned int
29#elif defined(INTEL_MKL)
30# define PREFIX
31# define Teuchos_fcd const char *
32#else /* Not CRAY_T3X or INTEL_CXML or INTEL_MKL */
33# define PREFIX
34# define Teuchos_fcd const char *
35#endif
36
37// Handle Complex types (we assume C++11 at a minumum
38// Microsoft C++11 apparently does not support float/double _Complex
39
40#if ( defined(_MSC_VER) )
41#define Teuchos_Complex_double_type_name std::complex<double>
42#define Teuchos_Complex_float_type_name std::complex<float>
43#else
44#define Teuchos_Complex_double_type_name double _Complex
45#define Teuchos_Complex_float_type_name float _Complex
46#endif
47
48
49/* B) Take care of of the link name case */
50
51
52#define DROTG_F77 F77_BLAS_MANGLE(drotg,DROTG)
53#define DROT_F77 F77_BLAS_MANGLE(drot,DROT)
54#define DASUM_F77 F77_BLAS_MANGLE(dasum,DASUM)
55#define DAXPY_F77 F77_BLAS_MANGLE(daxpy,DAXPY)
56#define DCOPY_F77 F77_BLAS_MANGLE(dcopy,DCOPY)
57#define DDOT_F77 F77_BLAS_MANGLE(ddot,DDOT)
58#define DNRM2_F77 F77_BLAS_MANGLE(dnrm2,DNRM2)
59#define DSCAL_F77 F77_BLAS_MANGLE(dscal,DSCAL)
60#define IDAMAX_F77 F77_BLAS_MANGLE(idamax,IDAMAX)
61#define DGEMV_F77 F77_BLAS_MANGLE(dgemv,DGEMV)
62#define DGER_F77 F77_BLAS_MANGLE(dger,DGER)
63#define DTRMV_F77 F77_BLAS_MANGLE(dtrmv,DTRMV)
64#define DGEMM_F77 F77_BLAS_MANGLE(dgemm,DGEMM)
65#define DSWAP_F77 F77_BLAS_MANGLE(dswap,DSWAP)
66#define DSYMM_F77 F77_BLAS_MANGLE(dsymm,DSYMM)
67#define DSYRK_F77 F77_BLAS_MANGLE(dsyrk,DSYRK)
68#define DTRMM_F77 F77_BLAS_MANGLE(dtrmm,DTRMM)
69#define DTRSM_F77 F77_BLAS_MANGLE(dtrsm,DTRSM)
70
71#ifdef HAVE_TEUCHOS_COMPLEX
72
73#define ZROTG_F77 F77_BLAS_MANGLE(zrotg,ZROTG)
74#define ZROT_F77 F77_BLAS_MANGLE(zrot,ZROT)
75#define ZASUM_F77 F77_BLAS_MANGLE(dzasum,DZASUM)
76#define ZAXPY_F77 F77_BLAS_MANGLE(zaxpy,ZAXPY)
77#define ZCOPY_F77 F77_BLAS_MANGLE(zcopy,ZCOPY)
78#define ZDOT_F77 F77_BLAS_MANGLE(zdotc,ZDOTC)
79#define ZNRM2_F77 F77_BLAS_MANGLE(dznrm2,DZNRM2)
80#define ZSCAL_F77 F77_BLAS_MANGLE(zscal,ZSCAL)
81#define IZAMAX_F77 F77_BLAS_MANGLE(izamax,IZAMAX)
82#define ZGEMV_F77 F77_BLAS_MANGLE(zgemv,ZGEMV)
83#define ZGER_F77 F77_BLAS_MANGLE(zgeru,ZGERU)
84#define ZTRMV_F77 F77_BLAS_MANGLE(ztrmv,ZTRMV)
85#define ZGEMM_F77 F77_BLAS_MANGLE(zgemm,ZGEMM)
86#define ZSWAP_F77 F77_BLAS_MANGLE(zswap,ZSWAP)
87#define ZSYMM_F77 F77_BLAS_MANGLE(zsymm,ZSYMM)
88#define ZSYRK_F77 F77_BLAS_MANGLE(zsyrk,ZSYRK)
89#define ZHERK_F77 F77_BLAS_MANGLE(zherk,ZHERK)
90#define ZTRMM_F77 F77_BLAS_MANGLE(ztrmm,ZTRMM)
91#define ZTRSM_F77 F77_BLAS_MANGLE(ztrsm,ZTRSM)
92
93#endif /* HAVE_TEUCHOS_COMPLEX */
94
95#define SROTG_F77 F77_BLAS_MANGLE(srotg,SROTG)
96#define SROT_F77 F77_BLAS_MANGLE(srot,SROT)
97#define SSCAL_F77 F77_BLAS_MANGLE(sscal,SSCAL)
98#define SCOPY_F77 F77_BLAS_MANGLE(scopy,SCOPY)
99#define SAXPY_F77 F77_BLAS_MANGLE(saxpy,SAXPY)
100#define SDOT_F77 F77_BLAS_MANGLE(sdot,SDOT)
101#define SNRM2_F77 F77_BLAS_MANGLE(snrm2,SNRM2)
102#define SASUM_F77 F77_BLAS_MANGLE(sasum,SASUM)
103#define ISAMAX_F77 F77_BLAS_MANGLE(isamax,ISAMAX)
104#define SGEMV_F77 F77_BLAS_MANGLE(sgemv,SGEMV)
105#define SGER_F77 F77_BLAS_MANGLE(sger,SGER)
106#define STRMV_F77 F77_BLAS_MANGLE(strmv,STRMV)
107#define SGEMM_F77 F77_BLAS_MANGLE(sgemm,SGEMM)
108#define SSWAP_F77 F77_BLAS_MANGLE(sswap,SSWAP)
109#define SSYMM_F77 F77_BLAS_MANGLE(ssymm,SSYMM)
110#define SSYRK_F77 F77_BLAS_MANGLE(ssyrk,SSYRK)
111#define STRMM_F77 F77_BLAS_MANGLE(strmm,STRMM)
112#define STRSM_F77 F77_BLAS_MANGLE(strsm,STRSM)
113
114#ifdef HAVE_TEUCHOS_COMPLEX
115
116#define CROTG_F77 F77_BLAS_MANGLE(crotg,CROTG)
117#define CROT_F77 F77_BLAS_MANGLE(crot,CROT)
118#define SCASUM_F77 F77_BLAS_MANGLE(scasum,SCASUM)
119#define CAXPY_F77 F77_BLAS_MANGLE(caxpy,CAXPY)
120#define CCOPY_F77 F77_BLAS_MANGLE(ccopy,CCOPY)
121#define CDOT_F77 F77_BLAS_MANGLE(cdotc,CDOTC)
122#define SCNRM2_F77 F77_BLAS_MANGLE(scnrm2,SCNRM2)
123#define CSCAL_F77 F77_BLAS_MANGLE(cscal,CSCAL)
124#define ICAMAX_F77 F77_BLAS_MANGLE(icamax,ICAMAX)
125#define CGEMV_F77 F77_BLAS_MANGLE(cgemv,CGEMV)
126#define CGER_F77 F77_BLAS_MANGLE(cgeru,CGERU)
127#define CTRMV_F77 F77_BLAS_MANGLE(ctrmv,CTRMV)
128#define CGEMM_F77 F77_BLAS_MANGLE(cgemm,CGEMM)
129#define CSWAP_F77 F77_BLAS_MANGLE(cswap,CSWAP)
130#define CSYMM_F77 F77_BLAS_MANGLE(csymm,CSYMM)
131#define CSYRK_F77 F77_BLAS_MANGLE(csyrk,CSYRK)
132#define CHERK_F77 F77_BLAS_MANGLE(cherk,CHERK)
133#define CTRMM_F77 F77_BLAS_MANGLE(ctrmm,CTRMM)
134#define CTRSM_F77 F77_BLAS_MANGLE(ctrsm,CTRSM)
135#define TEUCHOS_BLAS_CONVERT_COMPLEX_FORTRAN_TO_CXX(TYPE, Z) \
136 reinterpret_cast<std::complex<TYPE>&>(Z);
137// NOTE: The above is guaranteed to be okay given the C99 and C++11 standards
138
139#endif /* HAVE_TEUCHOS_COMPLEX */
140
141
142/* C) Define the function prototypes for all platforms! */
143
144#ifdef __cplusplus
145extern "C" {
146#endif
147
148
149/* Double precision BLAS 1 */
150void PREFIX DROTG_F77(double* da, double* db, double* c, double* s);
151void PREFIX DROT_F77(const int* n, double* dx, const int* incx, double* dy, const int* incy, double* c, double* s);
152double PREFIX DASUM_F77(const int* n, const double x[], const int* incx);
153void PREFIX DAXPY_F77(const int* n, const double* alpha, const double x[], const int* incx, double y[], const int* incy);
154void PREFIX DCOPY_F77(const int* n, const double *x, const int* incx, double *y, const int* incy);
155double PREFIX DDOT_F77(const int* n, const double x[], const int* incx, const double y[], const int* incy);
156double PREFIX DNRM2_F77(const int* n, const double x[], const int* incx);
157void PREFIX DSCAL_F77(const int* n, const double* alpha, double *x, const int* incx);
158void PREFIX DSWAP_F77(const int* const n, double* const x, const int* const incx,
159 double* const y, const int* const incy);
160int PREFIX IDAMAX_F77(const int* n, const double *x, const int* incx);
161
162/* Double std::complex precision BLAS 1 */
163#if defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus)
164
165# if defined(HAVE_COMPLEX_BLAS_PROBLEM)
166# if defined(HAVE_FIXABLE_COMPLEX_BLAS_PROBLEM)
167void PREFIX ZDOT_F77(std::complex<double> *ret, const int* n, const std::complex<double> x[], const int* incx, const std::complex<double> y[], const int* incy);
168# elif defined(HAVE_VECLIB_COMPLEX_BLAS)
169// no declarations; they're in cblas.h
170# include <vecLib/cblas.h>
171# else
172 // mfh 01 Feb 2013: If the code reaches this point, it means that
173 // some complex BLAS routines are broken, but there is no easy
174 // workaround. We deal with this in Teuchos_BLAS.cpp by
175 // reimplementing the offending routines.
176# endif // HAVE_COMPLEX_BLAS_PROBLEM
177# else // no problem
178Teuchos_Complex_double_type_name PREFIX ZDOT_F77(const int* n, const std::complex<double> x[], const int* incx, const std::complex<double> y[], const int* incy);
179# endif // defined(HAVE_COMPLEX_BLAS_PROBLEM)
180
181double PREFIX ZNRM2_F77(const int* n, const std::complex<double> x[], const int* incx);
182double PREFIX ZASUM_F77(const int* n, const std::complex<double> x[], const int* incx);
183void PREFIX ZROTG_F77(std::complex<double>* da, std::complex<double>* db, double* c, std::complex<double>* s);
184void PREFIX ZROT_F77(const int* n, std::complex<double>* dx, const int* incx, std::complex<double>* dy, const int* incy, double* c, std::complex<double>* s);
185void PREFIX ZAXPY_F77(const int* n, const std::complex<double>* alpha, const std::complex<double> x[], const int* incx, std::complex<double> y[], const int* incy);
186void PREFIX ZCOPY_F77(const int* n, const std::complex<double> *x, const int* incx, std::complex<double> *y, const int* incy);
187void PREFIX ZSCAL_F77(const int* n, const std::complex<double>* alpha, std::complex<double> *x, const int* incx);
188void PREFIX ZSWAP_F77(const int* const n, std::complex<double>* const x, const int* const incx,
189 std::complex<double>* const y, const int* const incy);
190int PREFIX IZAMAX_F77(const int* n, const std::complex<double> *x, const int* incx);
191
192#endif /* defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus) */
193
194/* Single precision BLAS 1 */
195#ifdef HAVE_TEUCHOS_BLASFLOAT
196# ifdef HAVE_TEUCHOS_BLASFLOAT_APPLE_VECLIB_BUGFIX
197# include <vecLib/cblas.h>
198# elif defined(HAVE_TEUCHOS_BLASFLOAT_DOUBLE_RETURN)
199double PREFIX SASUM_F77(const int* n, const float x[], const int* incx);
200double PREFIX SDOT_F77(const int* n, const float x[], const int* incx, const float y[], const int* incy);
201double PREFIX SNRM2_F77(const int* n, const float x[], const int* incx);
202# else
203float PREFIX SASUM_F77(const int* n, const float x[], const int* incx);
204float PREFIX SDOT_F77(const int* n, const float x[], const int* incx, const float y[], const int* incy);
205float PREFIX SNRM2_F77(const int* n, const float x[], const int* incx);
206# endif // which blasfloat
207#endif // ifdef blasfloat
208void PREFIX SROTG_F77(float* da, float* db, float* c, float* s);
209void PREFIX SROT_F77(const int* n, float* dx, const int* incx, float* dy, const int* incy, float* c, float* s);
210void PREFIX SAXPY_F77(const int* n, const float* alpha, const float x[], const int* incx, float y[], const int* incy);
211void PREFIX SCOPY_F77(const int* n, const float *x, const int* incx, float *y, const int* incy);
212void PREFIX SSCAL_F77(const int* n, const float* alpha, float *x, const int* incx);
213void PREFIX SSWAP_F77(const int* const n, float* const x, const int* const incx,
214 float* const y, const int* const incy);
215int PREFIX ISAMAX_F77(const int* n, const float *x, const int* incx);
216
217/* Single std::complex precision BLAS 1 */
218#if defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus)
219# if defined(HAVE_TEUCHOS_BLASFLOAT)
220# if defined(HAVE_TEUCHOS_BLASFLOAT_APPLE_VECLIB_BUGFIX)
221// no declarations; they're in cblas.h
222# include <vecLib/cblas.h>
223# elif defined(HAVE_TEUCHOS_BLASFLOAT_DOUBLE_RETURN)
224double PREFIX SCASUM_F77(const int* n, const std::complex<float> x[], const int* incx);
225double PREFIX SCNRM2_F77(const int* n, const std::complex<float> x[], const int* incx);
226# else
227float PREFIX SCASUM_F77(const int* n, const std::complex<float> x[], const int* incx);
228float PREFIX SCNRM2_F77(const int* n, const std::complex<float> x[], const int* incx);
229# endif // Whether or not we have the veclib bugfix
230#endif // defined(HAVE_TEUCHOS_BLASFLOAT)
231
232#if defined(HAVE_TEUCHOS_BLASFLOAT_APPLE_VECLIB_BUGFIX)
233// no declarations; they're in cblas.h
234#include <vecLib/cblas.h>
235#elif defined(HAVE_COMPLEX_BLAS_PROBLEM) && defined(HAVE_FIXABLE_COMPLEX_BLAS_PROBLEM)
236void PREFIX CDOT_F77(std::complex<float> *ret, const int* n, const std::complex<float> x[], const int* incx, const std::complex<float> y[], const int* incy);
237#elif defined(HAVE_TEUCHOS_BLASFLOAT)
238Teuchos_Complex_float_type_name PREFIX CDOT_F77(const int* n, const std::complex<float> x[], const int* incx, const std::complex<float> y[], const int* incy);
239#else
240// the code is literally in Teuchos_BLAS.cpp
241#endif
242
243void PREFIX CROTG_F77(std::complex<float>* da, std::complex<float>* db, float* c, std::complex<float>* s);
244void PREFIX CROT_F77(const int* n, std::complex<float>* dx, const int* incx, std::complex<float>* dy, const int* incy, float* c, std::complex<float>* s);
245void PREFIX CAXPY_F77(const int* n, const std::complex<float>* alpha, const std::complex<float> x[], const int* incx, std::complex<float> y[], const int* incy);
246void PREFIX CCOPY_F77(const int* n, const std::complex<float> *x, const int* incx, std::complex<float> *y, const int* incy);
247void PREFIX CSCAL_F77(const int* n, const std::complex<float>* alpha, std::complex<float> *x, const int* incx);
248void PREFIX CSWAP_F77(const int* const n, std::complex<float>* const x, const int* const incx,
249 std::complex<float>* const y, const int* const incy);
250int PREFIX ICAMAX_F77(const int* n, const std::complex<float> *x, const int* incx);
251
252#endif /* defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus) */
253
254/* Double precision BLAS 2 */
255void PREFIX DGEMV_F77(Teuchos_fcd, const int* m, const int* n, const double* alpha, const double A[], const int* lda,
256 const double x[], const int* incx, const double* beta, double y[], const int* incy);
257void PREFIX DTRMV_F77(Teuchos_fcd, Teuchos_fcd, Teuchos_fcd, const int *n,
258 const double *a, const int *lda, double *x, const int *incx);
259void PREFIX DGER_F77(const int *m, const int *n, const double *alpha, const double *x, const int *incx, const double *y,
260 const int *incy, double *a, const int *lda);
261
262/* Double precision BLAS 2 */
263#if defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus)
264
265void PREFIX ZGEMV_F77(Teuchos_fcd, const int* m, const int* n, const std::complex<double>* alpha, const std::complex<double> A[], const int* lda,
266 const std::complex<double> x[], const int* incx, const std::complex<double>* beta, std::complex<double> y[], const int* incy);
267void PREFIX ZTRMV_F77(Teuchos_fcd, Teuchos_fcd, Teuchos_fcd, const int *n,
268 const std::complex<double> *a, const int *lda, std::complex<double> *x, const int *incx);
269void PREFIX ZGER_F77(const int *m, const int *n, const std::complex<double> *alpha, const std::complex<double> *x, const int *incx, const std::complex<double> *y,
270 const int *incy, std::complex<double> *a, const int *lda);
271
272#endif /* defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus) */
273
274/* Single precision BLAS 2 */
275void PREFIX SGEMV_F77(Teuchos_fcd, const int* m, const int* n, const float* alpha, const float A[], const int* lda,
276 const float x[], const int* incx, const float* beta, float y[], const int* incy);
277void PREFIX STRMV_F77(Teuchos_fcd, Teuchos_fcd, Teuchos_fcd, const int *n,
278 const float *a, const int *lda, float *x, const int *incx);
279void PREFIX SGER_F77(const int *m, const int *n, const float *alpha, const float *x, const int *incx, const float *y,
280 const int *incy, float *a, const int *lda);
281
282/* Single std::complex precision BLAS 2 */
283#if defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus)
284
285void PREFIX CGEMV_F77(Teuchos_fcd, const int* m, const int* n, const std::complex<float>* alpha, const std::complex<float> A[], const int* lda,
286 const std::complex<float> x[], const int* incx, const std::complex<float>* beta, std::complex<float> y[], const int* incy);
287void PREFIX CTRMV_F77(Teuchos_fcd, Teuchos_fcd, Teuchos_fcd, const int *n,
288 const std::complex<float> *a, const int *lda, std::complex<float> *x, const int *incx);
289void PREFIX CGER_F77(const int *m, const int *n, const std::complex<float> *alpha, const std::complex<float> *x, const int *incx, const std::complex<float> *y,
290 const int *incy, std::complex<float> *a, const int *lda);
291
292#endif /* defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus) */
293
294/* Double precision BLAS 3 */
295void PREFIX DGEMM_F77(Teuchos_fcd, Teuchos_fcd, const int *m, const int *
296 n, const int *k, const double *alpha, const double *a, const int *lda,
297 const double *b, const int *ldb, const double *beta, double *c, const int *ldc);
298void PREFIX DSYMM_F77(Teuchos_fcd, Teuchos_fcd, const int *m, const int * n,
299 const double *alpha, const double *a, const int *lda,
300 const double *b, const int *ldb, const double *beta, double *c, const int *ldc);
301void PREFIX DSYRK_F77(Teuchos_fcd, Teuchos_fcd, const int *n, const int * k,
302 const double *alpha, const double *a, const int *lda,
303 const double *beta, double *c, const int *ldc);
304void PREFIX DTRMM_F77(Teuchos_fcd, Teuchos_fcd, Teuchos_fcd, Teuchos_fcd,
305 const int *m, const int *n, const double *alpha, const double *a, const int * lda, double *b, const int *ldb);
306void PREFIX DTRSM_F77(Teuchos_fcd, Teuchos_fcd, Teuchos_fcd, Teuchos_fcd,
307 const int *m, const int *n, const double *alpha, const double *a, const int *
308 lda, double *b, const int *ldb);
309
310/* Double std::complex precision BLAS 3 */
311#if defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus)
312
313void PREFIX ZGEMM_F77(Teuchos_fcd, Teuchos_fcd, const int *m, const int *
314 n, const int *k, const std::complex<double> *alpha, const std::complex<double> *a, const int *lda,
315 const std::complex<double> *b, const int *ldb, const std::complex<double> *beta, std::complex<double> *c, const int *ldc);
316void PREFIX ZSYMM_F77(Teuchos_fcd, Teuchos_fcd, const int *m, const int * n,
317 const std::complex<double> *alpha, const std::complex<double> *a, const int *lda,
318 const std::complex<double> *b, const int *ldb, const std::complex<double> *beta, std::complex<double> *c, const int *ldc);
319void PREFIX ZSYRK_F77(Teuchos_fcd, Teuchos_fcd, const int *n, const int * k,
320 const std::complex<double> *alpha, const std::complex<double> *a, const int *lda,
321 const std::complex<double> *beta, std::complex<double> *c, const int *ldc);
322void PREFIX ZHERK_F77(Teuchos_fcd, Teuchos_fcd, const int *n, const int * k,
323 const std::complex<double> *alpha, const std::complex<double> *a, const int *lda,
324 const std::complex<double> *beta, std::complex<double> *c, const int *ldc);
325void PREFIX ZTRMM_F77(Teuchos_fcd, Teuchos_fcd, Teuchos_fcd, Teuchos_fcd,
326 const int *m, const int *n, const std::complex<double> *alpha, const std::complex<double> *a, const int * lda, std::complex<double> *b, const int *ldb);
327void PREFIX ZTRSM_F77(Teuchos_fcd, Teuchos_fcd, Teuchos_fcd, Teuchos_fcd,
328 const int *m, const int *n, const std::complex<double> *alpha, const std::complex<double> *a, const int *
329 lda, std::complex<double> *b, const int *ldb);
330
331#endif /* defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus) */
332
333/* Single precision BLAS 3 */
334void PREFIX SGEMM_F77(Teuchos_fcd, Teuchos_fcd, const int *m, const int *
335 n, const int *k, const float *alpha, const float *a, const int *lda,
336 const float *b, const int *ldb, const float *beta, float *c, const int *ldc);
337void PREFIX SSYMM_F77(Teuchos_fcd, Teuchos_fcd, const int *m, const int * n,
338 const float *alpha, const float *a, const int *lda,
339 const float *b, const int *ldb, const float *beta, float *c, const int *ldc);
340void PREFIX SSYRK_F77(Teuchos_fcd, Teuchos_fcd, const int *n, const int * k,
341 const float *alpha, const float *a, const int *lda,
342 const float *beta, float *c, const int *ldc);
343void PREFIX STRMM_F77(Teuchos_fcd, Teuchos_fcd, Teuchos_fcd, Teuchos_fcd,
344 const int *m, const int *n, const float *alpha, const float *a, const int * lda, float *b, const int *ldb);
345void PREFIX STRSM_F77(Teuchos_fcd, Teuchos_fcd, Teuchos_fcd, Teuchos_fcd,
346 const int *m, const int *n, const float *alpha, const float *a, const int *
347 lda, float *b, const int *ldb);
348
349/* Single std::complex precision BLAS 3 */
350
351#if defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus)
352
353void PREFIX CGEMM_F77(Teuchos_fcd, Teuchos_fcd, const int *m, const int *
354 n, const int *k, const std::complex<float> *alpha, const std::complex<float> *a, const int *lda,
355 const std::complex<float> *b, const int *ldb, const std::complex<float> *beta, std::complex<float> *c, const int *ldc);
356void PREFIX CSYMM_F77(Teuchos_fcd, Teuchos_fcd, const int *m, const int * n,
357 const std::complex<float> *alpha, const std::complex<float> *a, const int *lda,
358 const std::complex<float> *b, const int *ldb, const std::complex<float> *beta, std::complex<float> *c, const int *ldc);
359void PREFIX CTRMM_F77(Teuchos_fcd, Teuchos_fcd, Teuchos_fcd, Teuchos_fcd,
360 const int *m, const int *n, const std::complex<float> *alpha, const std::complex<float> *a, const int * lda, std::complex<float> *b, const int *ldb);
361void PREFIX CSYRK_F77(Teuchos_fcd, Teuchos_fcd, const int *n, const int * k,
362 const std::complex<float> *alpha, const std::complex<float> *a, const int *lda,
363 const std::complex<float> *beta, std::complex<float> *c, const int *ldc);
364void PREFIX CHERK_F77(Teuchos_fcd, Teuchos_fcd, const int *n, const int * k,
365 const std::complex<float> *alpha, const std::complex<float> *a, const int *lda,
366 const std::complex<float> *beta, std::complex<float> *c, const int *ldc);
367void PREFIX CTRSM_F77(Teuchos_fcd, Teuchos_fcd, Teuchos_fcd, Teuchos_fcd,
368 const int *m, const int *n, const std::complex<float> *alpha, const std::complex<float> *a, const int *
369 lda, std::complex<float> *b, const int *ldb);
370
371#endif /* defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus) */
372
373#ifdef __cplusplus
374}
375#endif
376
377/* Don't leave a global macros called PREFIX or Teuchos_fcd laying around */
378
379#ifdef PREFIX
380#undef PREFIX
381#endif
382
383#ifdef Teuchos_fcd
384#undef Teuchos_fcd
385#endif
386
387#endif /* end of TEUCHOS_BLAS_WRAPPERS_HPP_ */
Teuchos header file which uses auto-configuration information to include necessary C++ headers.