Epetra Package Browser (Single Doxygen Collection)  Development
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Epetra_BLAS.cpp
Go to the documentation of this file.
1 
2 //@HEADER
3 // ************************************************************************
4 //
5 // Epetra: Linear Algebra Services Package
6 // Copyright 2011 Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
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 /* for INTEL_CXML, the second arg may need to be changed to 'one'. If so
44 the appropriate declaration of one will need to be added back into
45 functions that include the macro:
46 #if defined (INTEL_CXML)
47  unsigned int one=1;
48 #endif
49 */
50 
51 #ifdef CHAR_MACRO
52 #undef CHAR_MACRO
53 #endif
54 #if defined (INTEL_CXML)
55 #define CHAR_MACRO(char_var) &char_var, 1
56 #else
57 #define CHAR_MACRO(char_var) &char_var
58 #endif
59 
60 #include "Epetra_BLAS.h"
61 #include "Epetra_BLAS_wrappers.h"
62 
63 
64 //=============================================================================
65 float Epetra_BLAS::ASUM(const int N, const float * X, const int INCX) const {
66  return(SASUM_F77(&N, X, &INCX));
67 }
68 //=============================================================================
69 double Epetra_BLAS::ASUM(const int N, const double * X, const int INCX) const {
70  return(DASUM_F77(&N, X, &INCX));
71 }
72 //=============================================================================
73 float Epetra_BLAS::DOT(const int N, const float * X, const float * Y, const int INCX, const int INCY) const {
74  return(SDOT_F77(&N, X, &INCX, Y, &INCY));
75 }
76 //=============================================================================
77 double Epetra_BLAS::DOT(const int N, const double * X, const double * Y, const int INCX, const int INCY) const {
78  return(DDOT_F77(&N, X, &INCX, Y, &INCY));
79 }
80 //=============================================================================
81 float Epetra_BLAS::NRM2(const int N, const float * X, const int INCX) const {
82  return(SNRM2_F77(&N, X, &INCX));
83 }
84 //=============================================================================
85 double Epetra_BLAS::NRM2(const int N, const double * X, const int INCX) const {
86  return(DNRM2_F77(&N, X, &INCX));
87 }
88 //=============================================================================
89 void Epetra_BLAS::SCAL(const int N, const float ALPHA, float * X, const int INCX) const {
90  SSCAL_F77(&N, &ALPHA, X, &INCX);
91  return;
92 }
93 //=============================================================================
94 void Epetra_BLAS::SCAL(const int N, const double ALPHA, double * X, const int INCX) const {
95  DSCAL_F77(&N, &ALPHA, X, &INCX);
96  return;
97 }
98 //=============================================================================
99 void Epetra_BLAS::COPY(const int N, const float * X, float * Y, const int INCX, const int INCY) const {
100  SCOPY_F77(&N, X, &INCX, Y, &INCY);
101  return;
102 }
103 //=============================================================================
104 void Epetra_BLAS::COPY(const int N, const double * X, double * Y, const int INCX, const int INCY) const {
105  DCOPY_F77(&N, X, &INCX, Y, &INCY);
106  return;
107 }
108 //=============================================================================
109 int Epetra_BLAS::IAMAX(const int N, const float * X, const int INCX) const {
110  return(ISAMAX_F77(&N, X, &INCX)-1);// Note that we return base zero result
111 }
112 //=============================================================================
113 int Epetra_BLAS::IAMAX(const int N, const double * X, const int INCX) const {
114  return(IDAMAX_F77(&N, X, &INCX)-1);// Note that we return base zero result
115 }
116 //=============================================================================
117 void Epetra_BLAS::AXPY(const int N, const float ALPHA, const float * X, float * Y, const int INCX, const int INCY) const {
118  SAXPY_F77(&N, &ALPHA, X, &INCX, Y, &INCY);
119 }
120 //=============================================================================
121 void Epetra_BLAS::AXPY(const int N, const double ALPHA, const double * X, double * Y, const int INCX, const int INCY) const {
122  DAXPY_F77(&N, &ALPHA, X, &INCX, Y, &INCY);
123 }
124 //=============================================================================
125 void Epetra_BLAS::GEMV(const char TRANS, const int M, const int N,
126  const float ALPHA, const float * A, const int LDA, const float * X,
127  const float BETA, float * Y, const int INCX, const int INCY) const {
128  SGEMV_F77(CHAR_MACRO(TRANS), &M, &N, &ALPHA,
129  A, &LDA, X, &INCX, &BETA, Y, &INCY);
130 }
131 //=============================================================================
132 void Epetra_BLAS::GEMV(const char TRANS, const int M, const int N,
133  const double ALPHA, const double * A, const int LDA, const double * X,
134  const double BETA, double * Y, const int INCX, const int INCY) const {
135  DGEMV_F77(CHAR_MACRO(TRANS), &M, &N, &ALPHA,
136  A, &LDA, X, &INCX, &BETA, Y, &INCY);
137 }
138 
139 //=============================================================================
140 void Epetra_BLAS::GEMM(const char TRANSA, const char TRANSB, const int M, const int N, const int K,
141  const float ALPHA, const float * A, const int LDA, const float * B,
142  const int LDB, const float BETA, float * C, const int LDC) const {
143 
144  SGEMM_F77(CHAR_MACRO(TRANSA), CHAR_MACRO(TRANSB), &M, &N, &K, &ALPHA,
145  A, &LDA, B, &LDB, &BETA, C, &LDC);
146 }
147 
148 //=============================================================================
149 void Epetra_BLAS::GEMM(const char TRANSA, const char TRANSB, const int M, const int N, const int K,
150  const double ALPHA, const double * A, const int LDA, const double * B,
151  const int LDB, const double BETA, double * C, const int LDC) const {
152 
153  DGEMM_F77(CHAR_MACRO(TRANSA), CHAR_MACRO(TRANSB), &M, &N, &K, &ALPHA,
154  A, &LDA, B, &LDB, &BETA, C, &LDC);
155 }
156 //=============================================================================
157 void Epetra_BLAS::SYMM(const char SIDE, const char UPLO, const int M, const int N,
158  const float ALPHA, const float * A, const int LDA, const float * B,
159  const int LDB, const float BETA, float * C, const int LDC) const {
160 
161  SSYMM_F77(CHAR_MACRO(SIDE), CHAR_MACRO(UPLO), &M, &N, &ALPHA,
162  A, &LDA, B, &LDB, &BETA, C, &LDC);
163 }
164 
165 //=============================================================================
166 void Epetra_BLAS::SYMM(const char SIDE, const char UPLO, const int M, const int N,
167  const double ALPHA, const double * A, const int LDA, const double * B,
168  const int LDB, const double BETA, double * C, const int LDC) const {
169 
170  DSYMM_F77(CHAR_MACRO(SIDE), CHAR_MACRO(UPLO), &M, &N, &ALPHA,
171  A, &LDA, B, &LDB, &BETA, C, &LDC);
172 }
173 //=============================================================================
174 void Epetra_BLAS::TRMM(const char SIDE, const char UPLO, const char TRANSA, const char DIAG, const int M, const int N,
175  const float ALPHA, const float * A, const int LDA, float * B,
176  const int LDB) const {
177 
178  STRMM_F77(CHAR_MACRO(SIDE), CHAR_MACRO(UPLO), CHAR_MACRO(TRANSA), CHAR_MACRO(DIAG),
179  &M, &N, &ALPHA, A, &LDA, B, &LDB);
180 }
181 //=============================================================================
182 void Epetra_BLAS::TRMM(const char SIDE, const char UPLO, const char TRANSA, const char DIAG, const int M, const int N,
183  const double ALPHA, const double * A, const int LDA, double * B,
184  const int LDB) const {
185 
186  DTRMM_F77(CHAR_MACRO(SIDE), CHAR_MACRO(UPLO), CHAR_MACRO(TRANSA), CHAR_MACRO(DIAG),
187  &M, &N, &ALPHA, A, &LDA, B, &LDB);
188 }
189 //=============================================================================
190 void Epetra_BLAS::SYRK(const char UPLO, const char TRANS, const int N, const int K, const float ALPHA, const float *A,
191  const int LDA, const float BETA, float *C, const int LDC) const{
192  SSYRK_F77(CHAR_MACRO(UPLO), CHAR_MACRO(TRANS), &N, &K, &ALPHA, A, &LDA, &BETA, C, &LDC);
193 }
194 //=============================================================================
195 void Epetra_BLAS::SYRK(const char UPLO, const char TRANS, const int N, const int K, const double ALPHA, const double *A,
196  const int LDA, const double BETA, double *C, const int LDC) const{
197  DSYRK_F77(CHAR_MACRO(UPLO), CHAR_MACRO(TRANS), &N, &K, &ALPHA, A, &LDA, &BETA, C, &LDC);
198 }
#define SGEMM_F77
float PREFIX SDOT_F77(const int *n, const float x[], const int *incx, const float y[], const int *incy)
void GEMM(const char TRANSA, const char TRANSB, const int M, const int N, const int K, const float ALPHA, const float *A, const int LDA, const float *B, const int LDB, const float BETA, float *C, const int LDC) const
Epetra_BLAS matrix-matrix multiply function (SGEMM)
#define DCOPY_F77
#define SSYRK_F77
#define CHAR_MACRO(char_var)
Definition: Epetra_BLAS.cpp:57
void AXPY(const int N, const float ALPHA, const float *X, float *Y, const int INCX=1, const int INCY=1) const
Epetra_BLAS vector update function (SAXPY)
#define ISAMAX_F77
float NRM2(const int N, const float *X, const int INCX=1) const
Epetra_BLAS norm function (SNRM2).
Definition: Epetra_BLAS.cpp:81
void TRMM(const char SIDE, const char UPLO, const char TRANSA, const char DIAG, const int M, const int N, const float ALPHA, const float *A, const int LDA, float *B, const int LDB) const
Epetra_BLAS triangular matrix-matrix multiply function (STRMM)
#define STRMM_F77
void GEMV(const char TRANS, const int M, const int N, const float ALPHA, const float *A, const int LDA, const float *X, const float BETA, float *Y, const int INCX=1, const int INCY=1) const
Epetra_BLAS matrix-vector multiply function (SGEMV)
#define DAXPY_F77
#define SSYMM_F77
#define DDOT_F77
#define SCOPY_F77
#define DTRMM_F77
#define DASUM_F77
#define DNRM2_F77
#define SSCAL_F77
float PREFIX SNRM2_F77(const int *n, const float x[], const int *incx)
void SYRK(const char UPLO, const char TRANS, const int N, const int K, const float ALPHA, const float *A, const int LDA, const float BETA, float *C, const int LDC) const
Eperta_BLAS symetric rank k funtion (ssyrk)
float PREFIX SASUM_F77(const int *n, const float x[], const int *incx)
#define DGEMV_F77
float ASUM(const int N, const float *X, const int INCX=1) const
Epetra_BLAS one norm function (SASUM).
Definition: Epetra_BLAS.cpp:65
#define DSCAL_F77
void SYMM(const char SIDE, const char UPLO, const int M, const int N, const float ALPHA, const float *A, const int LDA, const float *B, const int LDB, const float BETA, float *C, const int LDC) const
Epetra_BLAS symmetric matrix-matrix multiply function (SSYMM)
#define DGEMM_F77
void SCAL(const int N, const float ALPHA, float *X, const int INCX=1) const
Epetra_BLAS vector scale function (SSCAL)
Definition: Epetra_BLAS.cpp:89
int IAMAX(const int N, const float *X, const int INCX=1) const
Epetra_BLAS arg maximum of absolute value function (ISAMAX)
#define IDAMAX_F77
float DOT(const int N, const float *X, const float *Y, const int INCX=1, const int INCY=1) const
Epetra_BLAS dot product function (SDOT).
Definition: Epetra_BLAS.cpp:73
void COPY(const int N, const float *X, float *Y, const int INCX=1, const int INCY=1) const
Epetra_BLAS vector copy function (SCOPY)
Definition: Epetra_BLAS.cpp:99
#define SGEMV_F77
#define DSYRK_F77
#define DSYMM_F77
#define SAXPY_F77