Teuchos Package Browser (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Teuchos_BLAS_wrappers.hpp
Go to the documentation of this file.
1 /*
2 // @HEADER
3 // ***********************************************************************
4 //
5 // Teuchos: Common Tools Package
6 // Copyright (2004) Sandia Corporation
7 //
8 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
9 // license for use of this work by or on behalf of the U.S. Government.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
39 //
40 // ***********************************************************************
41 // @HEADER
42 */
43 
44 #ifndef _TEUCHOS_BLAS_WRAPPERS_HPP_
45 #define _TEUCHOS_BLAS_WRAPPERS_HPP_
46 
47 #include "Teuchos_ConfigDefs.hpp"
48 #ifdef _MSC_VER
49 /* disable warning for C-linkage returning complex class */
50 #pragma warning ( disable : 4190 )
51 #endif
52 
58 /* A) Define PREFIX and Teuchos_fcd based on platform. */
59 
60 #if defined(INTEL_CXML)
61 # define PREFIX __stdcall
62 # define Teuchos_fcd const char *, unsigned int
63 #elif defined(INTEL_MKL)
64 # define PREFIX
65 # define Teuchos_fcd const char *
66 #else /* Not CRAY_T3X or INTEL_CXML or INTEL_MKL */
67 # define PREFIX
68 # define Teuchos_fcd const char *
69 #endif
70 
71 // Handle Complex types (we assume C++11 at a minumum
72 // Microsoft C++11 apparently does not support float/double _Complex
73 
74 #if ( defined(_MSC_VER) )
75 #define Teuchos_Complex_double_type_name std::complex<double>
76 #define Teuchos_Complex_float_type_name std::complex<float>
77 #else
78 #define Teuchos_Complex_double_type_name double _Complex
79 #define Teuchos_Complex_float_type_name float _Complex
80 #endif
81 
82 
83 /* B) Take care of of the link name case */
84 
85 
86 #define DROTG_F77 F77_BLAS_MANGLE(drotg,DROTG)
87 #define DROT_F77 F77_BLAS_MANGLE(drot,DROT)
88 #define DASUM_F77 F77_BLAS_MANGLE(dasum,DASUM)
89 #define DAXPY_F77 F77_BLAS_MANGLE(daxpy,DAXPY)
90 #define DCOPY_F77 F77_BLAS_MANGLE(dcopy,DCOPY)
91 #define DDOT_F77 F77_BLAS_MANGLE(ddot,DDOT)
92 #define DNRM2_F77 F77_BLAS_MANGLE(dnrm2,DNRM2)
93 #define DSCAL_F77 F77_BLAS_MANGLE(dscal,DSCAL)
94 #define IDAMAX_F77 F77_BLAS_MANGLE(idamax,IDAMAX)
95 #define DGEMV_F77 F77_BLAS_MANGLE(dgemv,DGEMV)
96 #define DGER_F77 F77_BLAS_MANGLE(dger,DGER)
97 #define DTRMV_F77 F77_BLAS_MANGLE(dtrmv,DTRMV)
98 #define DGEMM_F77 F77_BLAS_MANGLE(dgemm,DGEMM)
99 #define DSWAP_F77 F77_BLAS_MANGLE(dswap,DSWAP)
100 #define DSYMM_F77 F77_BLAS_MANGLE(dsymm,DSYMM)
101 #define DSYRK_F77 F77_BLAS_MANGLE(dsyrk,DSYRK)
102 #define DTRMM_F77 F77_BLAS_MANGLE(dtrmm,DTRMM)
103 #define DTRSM_F77 F77_BLAS_MANGLE(dtrsm,DTRSM)
104 
105 #ifdef HAVE_TEUCHOS_COMPLEX
106 
107 #define ZROTG_F77 F77_BLAS_MANGLE(zrotg,ZROTG)
108 #define ZROT_F77 F77_BLAS_MANGLE(zrot,ZROT)
109 #define ZASUM_F77 F77_BLAS_MANGLE(dzasum,DZASUM)
110 #define ZAXPY_F77 F77_BLAS_MANGLE(zaxpy,ZAXPY)
111 #define ZCOPY_F77 F77_BLAS_MANGLE(zcopy,ZCOPY)
112 #define ZDOT_F77 F77_BLAS_MANGLE(zdotc,ZDOTC)
113 #define ZNRM2_F77 F77_BLAS_MANGLE(dznrm2,DZNRM2)
114 #define ZSCAL_F77 F77_BLAS_MANGLE(zscal,ZSCAL)
115 #define IZAMAX_F77 F77_BLAS_MANGLE(izamax,IZAMAX)
116 #define ZGEMV_F77 F77_BLAS_MANGLE(zgemv,ZGEMV)
117 #define ZGER_F77 F77_BLAS_MANGLE(zgeru,ZGERU)
118 #define ZTRMV_F77 F77_BLAS_MANGLE(ztrmv,ZTRMV)
119 #define ZGEMM_F77 F77_BLAS_MANGLE(zgemm,ZGEMM)
120 #define ZSWAP_F77 F77_BLAS_MANGLE(zswap,ZSWAP)
121 #define ZSYMM_F77 F77_BLAS_MANGLE(zsymm,ZSYMM)
122 #define ZSYRK_F77 F77_BLAS_MANGLE(zsyrk,ZSYRK)
123 #define ZHERK_F77 F77_BLAS_MANGLE(zherk,ZHERK)
124 #define ZTRMM_F77 F77_BLAS_MANGLE(ztrmm,ZTRMM)
125 #define ZTRSM_F77 F77_BLAS_MANGLE(ztrsm,ZTRSM)
126 
127 #endif /* HAVE_TEUCHOS_COMPLEX */
128 
129 #define SROTG_F77 F77_BLAS_MANGLE(srotg,SROTG)
130 #define SROT_F77 F77_BLAS_MANGLE(srot,SROT)
131 #define SSCAL_F77 F77_BLAS_MANGLE(sscal,SSCAL)
132 #define SCOPY_F77 F77_BLAS_MANGLE(scopy,SCOPY)
133 #define SAXPY_F77 F77_BLAS_MANGLE(saxpy,SAXPY)
134 #define SDOT_F77 F77_BLAS_MANGLE(sdot,SDOT)
135 #define SNRM2_F77 F77_BLAS_MANGLE(snrm2,SNRM2)
136 #define SASUM_F77 F77_BLAS_MANGLE(sasum,SASUM)
137 #define ISAMAX_F77 F77_BLAS_MANGLE(isamax,ISAMAX)
138 #define SGEMV_F77 F77_BLAS_MANGLE(sgemv,SGEMV)
139 #define SGER_F77 F77_BLAS_MANGLE(sger,SGER)
140 #define STRMV_F77 F77_BLAS_MANGLE(strmv,STRMV)
141 #define SGEMM_F77 F77_BLAS_MANGLE(sgemm,SGEMM)
142 #define SSWAP_F77 F77_BLAS_MANGLE(sswap,SSWAP)
143 #define SSYMM_F77 F77_BLAS_MANGLE(ssymm,SSYMM)
144 #define SSYRK_F77 F77_BLAS_MANGLE(ssyrk,SSYRK)
145 #define STRMM_F77 F77_BLAS_MANGLE(strmm,STRMM)
146 #define STRSM_F77 F77_BLAS_MANGLE(strsm,STRSM)
147 
148 #ifdef HAVE_TEUCHOS_COMPLEX
149 
150 #define CROTG_F77 F77_BLAS_MANGLE(crotg,CROTG)
151 #define CROT_F77 F77_BLAS_MANGLE(crot,CROT)
152 #define SCASUM_F77 F77_BLAS_MANGLE(scasum,SCASUM)
153 #define CAXPY_F77 F77_BLAS_MANGLE(caxpy,CAXPY)
154 #define CCOPY_F77 F77_BLAS_MANGLE(ccopy,CCOPY)
155 #define CDOT_F77 F77_BLAS_MANGLE(cdotc,CDOTC)
156 #define SCNRM2_F77 F77_BLAS_MANGLE(scnrm2,SCNRM2)
157 #define CSCAL_F77 F77_BLAS_MANGLE(cscal,CSCAL)
158 #define ICAMAX_F77 F77_BLAS_MANGLE(icamax,ICAMAX)
159 #define CGEMV_F77 F77_BLAS_MANGLE(cgemv,CGEMV)
160 #define CGER_F77 F77_BLAS_MANGLE(cgeru,CGERU)
161 #define CTRMV_F77 F77_BLAS_MANGLE(ctrmv,CTRMV)
162 #define CGEMM_F77 F77_BLAS_MANGLE(cgemm,CGEMM)
163 #define CSWAP_F77 F77_BLAS_MANGLE(cswap,CSWAP)
164 #define CSYMM_F77 F77_BLAS_MANGLE(csymm,CSYMM)
165 #define CSYRK_F77 F77_BLAS_MANGLE(csyrk,CSYRK)
166 #define CHERK_F77 F77_BLAS_MANGLE(cherk,CHERK)
167 #define CTRMM_F77 F77_BLAS_MANGLE(ctrmm,CTRMM)
168 #define CTRSM_F77 F77_BLAS_MANGLE(ctrsm,CTRSM)
169 #define TEUCHOS_BLAS_CONVERT_COMPLEX_FORTRAN_TO_CXX(TYPE, Z) \
170  reinterpret_cast<std::complex<TYPE>&>(Z);
171 // NOTE: The above is guaranteed to be okay given the C99 and C++11 standards
172 
173 #endif /* HAVE_TEUCHOS_COMPLEX */
174 
175 
176 /* C) Define the function prototypes for all platforms! */
177 
178 #ifdef __cplusplus
179 extern "C" {
180 #endif
181 
182 
183 /* Double precision BLAS 1 */
184 void PREFIX DROTG_F77(double* da, double* db, double* c, double* s);
185 void PREFIX DROT_F77(const int* n, double* dx, const int* incx, double* dy, const int* incy, double* c, double* s);
186 double PREFIX DASUM_F77(const int* n, const double x[], const int* incx);
187 void PREFIX DAXPY_F77(const int* n, const double* alpha, const double x[], const int* incx, double y[], const int* incy);
188 void PREFIX DCOPY_F77(const int* n, const double *x, const int* incx, double *y, const int* incy);
189 double PREFIX DDOT_F77(const int* n, const double x[], const int* incx, const double y[], const int* incy);
190 double PREFIX DNRM2_F77(const int* n, const double x[], const int* incx);
191 void PREFIX DSCAL_F77(const int* n, const double* alpha, double *x, const int* incx);
192 void PREFIX DSWAP_F77(const int* const n, double* const x, const int* const incx,
193  double* const y, const int* const incy);
194 int PREFIX IDAMAX_F77(const int* n, const double *x, const int* incx);
195 
196 /* Double std::complex precision BLAS 1 */
197 #if defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus)
198 
199 # if defined(HAVE_COMPLEX_BLAS_PROBLEM)
200 # if defined(HAVE_FIXABLE_COMPLEX_BLAS_PROBLEM)
201 void 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);
202 # elif defined(HAVE_VECLIB_COMPLEX_BLAS)
203 // no declarations; they're in cblas.h
204 # include <vecLib/cblas.h>
205 # else
206  // mfh 01 Feb 2013: If the code reaches this point, it means that
207  // some complex BLAS routines are broken, but there is no easy
208  // workaround. We deal with this in Teuchos_BLAS.cpp by
209  // reimplementing the offending routines.
210 # endif // HAVE_COMPLEX_BLAS_PROBLEM
211 # else // no problem
212 Teuchos_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);
213 # endif // defined(HAVE_COMPLEX_BLAS_PROBLEM)
214 
215 double PREFIX ZNRM2_F77(const int* n, const std::complex<double> x[], const int* incx);
216 double PREFIX ZASUM_F77(const int* n, const std::complex<double> x[], const int* incx);
217 void PREFIX ZROTG_F77(std::complex<double>* da, std::complex<double>* db, double* c, std::complex<double>* s);
218 void 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);
219 void 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);
220 void PREFIX ZCOPY_F77(const int* n, const std::complex<double> *x, const int* incx, std::complex<double> *y, const int* incy);
221 void PREFIX ZSCAL_F77(const int* n, const std::complex<double>* alpha, std::complex<double> *x, const int* incx);
222 void PREFIX ZSWAP_F77(const int* const n, std::complex<double>* const x, const int* const incx,
223  std::complex<double>* const y, const int* const incy);
224 int PREFIX IZAMAX_F77(const int* n, const std::complex<double> *x, const int* incx);
225 
226 #endif /* defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus) */
227 
228 /* Single precision BLAS 1 */
229 #ifdef HAVE_TEUCHOS_BLASFLOAT
230 # ifdef HAVE_TEUCHOS_BLASFLOAT_APPLE_VECLIB_BUGFIX
231 # include <vecLib/cblas.h>
232 # elif defined(HAVE_TEUCHOS_BLASFLOAT_DOUBLE_RETURN)
233 double PREFIX SASUM_F77(const int* n, const float x[], const int* incx);
234 double PREFIX SDOT_F77(const int* n, const float x[], const int* incx, const float y[], const int* incy);
235 double PREFIX SNRM2_F77(const int* n, const float x[], const int* incx);
236 # else
237 float PREFIX SASUM_F77(const int* n, const float x[], const int* incx);
238 float PREFIX SDOT_F77(const int* n, const float x[], const int* incx, const float y[], const int* incy);
239 float PREFIX SNRM2_F77(const int* n, const float x[], const int* incx);
240 # endif // which blasfloat
241 #endif // ifdef blasfloat
242 void PREFIX SROTG_F77(float* da, float* db, float* c, float* s);
243 void PREFIX SROT_F77(const int* n, float* dx, const int* incx, float* dy, const int* incy, float* c, float* s);
244 void PREFIX SAXPY_F77(const int* n, const float* alpha, const float x[], const int* incx, float y[], const int* incy);
245 void PREFIX SCOPY_F77(const int* n, const float *x, const int* incx, float *y, const int* incy);
246 void PREFIX SSCAL_F77(const int* n, const float* alpha, float *x, const int* incx);
247 void PREFIX SSWAP_F77(const int* const n, float* const x, const int* const incx,
248  float* const y, const int* const incy);
249 int PREFIX ISAMAX_F77(const int* n, const float *x, const int* incx);
250 
251 /* Single std::complex precision BLAS 1 */
252 #if defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus)
253 # if defined(HAVE_TEUCHOS_BLASFLOAT)
254 # if defined(HAVE_TEUCHOS_BLASFLOAT_APPLE_VECLIB_BUGFIX)
255 // no declarations; they're in cblas.h
256 # include <vecLib/cblas.h>
257 # elif defined(HAVE_TEUCHOS_BLASFLOAT_DOUBLE_RETURN)
258 double PREFIX SCASUM_F77(const int* n, const std::complex<float> x[], const int* incx);
259 double PREFIX SCNRM2_F77(const int* n, const std::complex<float> x[], const int* incx);
260 # else
261 float PREFIX SCASUM_F77(const int* n, const std::complex<float> x[], const int* incx);
262 float PREFIX SCNRM2_F77(const int* n, const std::complex<float> x[], const int* incx);
263 # endif // Whether or not we have the veclib bugfix
264 #endif // defined(HAVE_TEUCHOS_BLASFLOAT)
265 
266 #if defined(HAVE_TEUCHOS_BLASFLOAT_APPLE_VECLIB_BUGFIX)
267 // no declarations; they're in cblas.h
268 #include <vecLib/cblas.h>
269 #elif defined(HAVE_COMPLEX_BLAS_PROBLEM) && defined(HAVE_FIXABLE_COMPLEX_BLAS_PROBLEM)
270 void 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);
271 #elif defined(HAVE_TEUCHOS_BLASFLOAT)
272 Teuchos_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);
273 #else
274 // the code is literally in Teuchos_BLAS.cpp
275 #endif
276 
277 void PREFIX CROTG_F77(std::complex<float>* da, std::complex<float>* db, float* c, std::complex<float>* s);
278 void 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);
279 void 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);
280 void PREFIX CCOPY_F77(const int* n, const std::complex<float> *x, const int* incx, std::complex<float> *y, const int* incy);
281 void PREFIX CSCAL_F77(const int* n, const std::complex<float>* alpha, std::complex<float> *x, const int* incx);
282 void PREFIX CSWAP_F77(const int* const n, std::complex<float>* const x, const int* const incx,
283  std::complex<float>* const y, const int* const incy);
284 int PREFIX ICAMAX_F77(const int* n, const std::complex<float> *x, const int* incx);
285 
286 #endif /* defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus) */
287 
288 /* Double precision BLAS 2 */
289 void PREFIX DGEMV_F77(Teuchos_fcd, const int* m, const int* n, const double* alpha, const double A[], const int* lda,
290  const double x[], const int* incx, const double* beta, double y[], const int* incy);
291 void PREFIX DTRMV_F77(Teuchos_fcd, Teuchos_fcd, Teuchos_fcd, const int *n,
292  const double *a, const int *lda, double *x, const int *incx);
293 void PREFIX DGER_F77(const int *m, const int *n, const double *alpha, const double *x, const int *incx, const double *y,
294  const int *incy, double *a, const int *lda);
295 
296 /* Double precision BLAS 2 */
297 #if defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus)
298 
299 void PREFIX ZGEMV_F77(Teuchos_fcd, const int* m, const int* n, const std::complex<double>* alpha, const std::complex<double> A[], const int* lda,
300  const std::complex<double> x[], const int* incx, const std::complex<double>* beta, std::complex<double> y[], const int* incy);
301 void PREFIX ZTRMV_F77(Teuchos_fcd, Teuchos_fcd, Teuchos_fcd, const int *n,
302  const std::complex<double> *a, const int *lda, std::complex<double> *x, const int *incx);
303 void 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,
304  const int *incy, std::complex<double> *a, const int *lda);
305 
306 #endif /* defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus) */
307 
308 /* Single precision BLAS 2 */
309 void PREFIX SGEMV_F77(Teuchos_fcd, const int* m, const int* n, const float* alpha, const float A[], const int* lda,
310  const float x[], const int* incx, const float* beta, float y[], const int* incy);
311 void PREFIX STRMV_F77(Teuchos_fcd, Teuchos_fcd, Teuchos_fcd, const int *n,
312  const float *a, const int *lda, float *x, const int *incx);
313 void PREFIX SGER_F77(const int *m, const int *n, const float *alpha, const float *x, const int *incx, const float *y,
314  const int *incy, float *a, const int *lda);
315 
316 /* Single std::complex precision BLAS 2 */
317 #if defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus)
318 
319 void PREFIX CGEMV_F77(Teuchos_fcd, const int* m, const int* n, const std::complex<float>* alpha, const std::complex<float> A[], const int* lda,
320  const std::complex<float> x[], const int* incx, const std::complex<float>* beta, std::complex<float> y[], const int* incy);
321 void PREFIX CTRMV_F77(Teuchos_fcd, Teuchos_fcd, Teuchos_fcd, const int *n,
322  const std::complex<float> *a, const int *lda, std::complex<float> *x, const int *incx);
323 void 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,
324  const int *incy, std::complex<float> *a, const int *lda);
325 
326 #endif /* defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus) */
327 
328 /* Double precision BLAS 3 */
329 void PREFIX DGEMM_F77(Teuchos_fcd, Teuchos_fcd, const int *m, const int *
330  n, const int *k, const double *alpha, const double *a, const int *lda,
331  const double *b, const int *ldb, const double *beta, double *c, const int *ldc);
332 void PREFIX DSYMM_F77(Teuchos_fcd, Teuchos_fcd, const int *m, const int * n,
333  const double *alpha, const double *a, const int *lda,
334  const double *b, const int *ldb, const double *beta, double *c, const int *ldc);
335 void PREFIX DSYRK_F77(Teuchos_fcd, Teuchos_fcd, const int *n, const int * k,
336  const double *alpha, const double *a, const int *lda,
337  const double *beta, double *c, const int *ldc);
339  const int *m, const int *n, const double *alpha, const double *a, const int * lda, double *b, const int *ldb);
341  const int *m, const int *n, const double *alpha, const double *a, const int *
342  lda, double *b, const int *ldb);
343 
344 /* Double std::complex precision BLAS 3 */
345 #if defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus)
346 
347 void PREFIX ZGEMM_F77(Teuchos_fcd, Teuchos_fcd, const int *m, const int *
348  n, const int *k, const std::complex<double> *alpha, const std::complex<double> *a, const int *lda,
349  const std::complex<double> *b, const int *ldb, const std::complex<double> *beta, std::complex<double> *c, const int *ldc);
350 void PREFIX ZSYMM_F77(Teuchos_fcd, Teuchos_fcd, const int *m, const int * n,
351  const std::complex<double> *alpha, const std::complex<double> *a, const int *lda,
352  const std::complex<double> *b, const int *ldb, const std::complex<double> *beta, std::complex<double> *c, const int *ldc);
353 void PREFIX ZSYRK_F77(Teuchos_fcd, Teuchos_fcd, const int *n, const int * k,
354  const std::complex<double> *alpha, const std::complex<double> *a, const int *lda,
355  const std::complex<double> *beta, std::complex<double> *c, const int *ldc);
356 void PREFIX ZHERK_F77(Teuchos_fcd, Teuchos_fcd, const int *n, const int * k,
357  const std::complex<double> *alpha, const std::complex<double> *a, const int *lda,
358  const std::complex<double> *beta, std::complex<double> *c, const int *ldc);
360  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);
362  const int *m, const int *n, const std::complex<double> *alpha, const std::complex<double> *a, const int *
363  lda, std::complex<double> *b, const int *ldb);
364 
365 #endif /* defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus) */
366 
367 /* Single precision BLAS 3 */
368 void PREFIX SGEMM_F77(Teuchos_fcd, Teuchos_fcd, const int *m, const int *
369  n, const int *k, const float *alpha, const float *a, const int *lda,
370  const float *b, const int *ldb, const float *beta, float *c, const int *ldc);
371 void PREFIX SSYMM_F77(Teuchos_fcd, Teuchos_fcd, const int *m, const int * n,
372  const float *alpha, const float *a, const int *lda,
373  const float *b, const int *ldb, const float *beta, float *c, const int *ldc);
374 void PREFIX SSYRK_F77(Teuchos_fcd, Teuchos_fcd, const int *n, const int * k,
375  const float *alpha, const float *a, const int *lda,
376  const float *beta, float *c, const int *ldc);
378  const int *m, const int *n, const float *alpha, const float *a, const int * lda, float *b, const int *ldb);
380  const int *m, const int *n, const float *alpha, const float *a, const int *
381  lda, float *b, const int *ldb);
382 
383 /* Single std::complex precision BLAS 3 */
384 
385 #if defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus)
386 
387 void PREFIX CGEMM_F77(Teuchos_fcd, Teuchos_fcd, const int *m, const int *
388  n, const int *k, const std::complex<float> *alpha, const std::complex<float> *a, const int *lda,
389  const std::complex<float> *b, const int *ldb, const std::complex<float> *beta, std::complex<float> *c, const int *ldc);
390 void PREFIX CSYMM_F77(Teuchos_fcd, Teuchos_fcd, const int *m, const int * n,
391  const std::complex<float> *alpha, const std::complex<float> *a, const int *lda,
392  const std::complex<float> *b, const int *ldb, const std::complex<float> *beta, std::complex<float> *c, const int *ldc);
394  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);
395 void PREFIX CSYRK_F77(Teuchos_fcd, Teuchos_fcd, const int *n, const int * k,
396  const std::complex<float> *alpha, const std::complex<float> *a, const int *lda,
397  const std::complex<float> *beta, std::complex<float> *c, const int *ldc);
398 void PREFIX CHERK_F77(Teuchos_fcd, Teuchos_fcd, const int *n, const int * k,
399  const std::complex<float> *alpha, const std::complex<float> *a, const int *lda,
400  const std::complex<float> *beta, std::complex<float> *c, const int *ldc);
402  const int *m, const int *n, const std::complex<float> *alpha, const std::complex<float> *a, const int *
403  lda, std::complex<float> *b, const int *ldb);
404 
405 #endif /* defined(HAVE_TEUCHOS_COMPLEX) && defined(__cplusplus) */
406 
407 #ifdef __cplusplus
408 }
409 #endif
410 
411 /* Don't leave a global macros called PREFIX or Teuchos_fcd laying around */
412 
413 #ifdef PREFIX
414 #undef PREFIX
415 #endif
416 
417 #ifdef Teuchos_fcd
418 #undef Teuchos_fcd
419 #endif
420 
421 #endif /* end of TEUCHOS_BLAS_WRAPPERS_HPP_ */
#define Teuchos_Complex_double_type_name
#define PREFIX
#define STRMM_F77
#define SSCAL_F77
#define DTRMV_F77
#define IDAMAX_F77
#define Teuchos_Complex_float_type_name
#define DROTG_F77
Teuchos header file which uses auto-configuration information to include necessary C++ headers...
#define SROTG_F77
#define SGER_F77
#define STRSM_F77
#define DAXPY_F77
#define SDOT_F77
#define SROT_F77
#define Teuchos_fcd
#define SAXPY_F77
#define DSYMM_F77
#define DGER_F77
#define SSYRK_F77
#define DNRM2_F77
#define DSWAP_F77
#define DTRMM_F77
#define SGEMM_F77
#define DGEMV_F77
#define SASUM_F77
#define SCOPY_F77
#define DCOPY_F77
#define ISAMAX_F77
#define DSCAL_F77
#define DSYRK_F77
#define SSYMM_F77
#define DTRSM_F77
#define SSWAP_F77
#define DASUM_F77
#define DDOT_F77
#define SGEMV_F77
#define DGEMM_F77
#define DROT_F77
#define STRMV_F77
#define SNRM2_F77