Chris@49
|
1 // Copyright (C) 2008-2012 NICTA (www.nicta.com.au)
|
Chris@49
|
2 // Copyright (C) 2008-2012 Conrad Sanderson
|
Chris@49
|
3 //
|
Chris@49
|
4 // This Source Code Form is subject to the terms of the Mozilla Public
|
Chris@49
|
5 // License, v. 2.0. If a copy of the MPL was not distributed with this
|
Chris@49
|
6 // file, You can obtain one at http://mozilla.org/MPL/2.0/.
|
Chris@49
|
7
|
Chris@49
|
8
|
Chris@49
|
9 #ifdef ARMA_USE_ATLAS
|
Chris@49
|
10
|
Chris@49
|
11
|
Chris@49
|
12 //! \namespace atlas namespace for ATLAS functions (imported from the global namespace)
|
Chris@49
|
13 namespace atlas
|
Chris@49
|
14 {
|
Chris@49
|
15
|
Chris@49
|
16 template<typename eT>
|
Chris@49
|
17 inline static const eT& tmp_real(const eT& X) { return X; }
|
Chris@49
|
18
|
Chris@49
|
19 template<typename T>
|
Chris@49
|
20 inline static const T tmp_real(const std::complex<T>& X) { return X.real(); }
|
Chris@49
|
21
|
Chris@49
|
22
|
Chris@49
|
23
|
Chris@49
|
24 template<typename eT>
|
Chris@49
|
25 arma_inline
|
Chris@49
|
26 eT
|
Chris@49
|
27 cblas_dot(const int N, const eT* X, const eT* Y)
|
Chris@49
|
28 {
|
Chris@49
|
29 arma_type_check((is_supported_blas_type<eT>::value == false));
|
Chris@49
|
30
|
Chris@49
|
31 if(is_float<eT>::value == true)
|
Chris@49
|
32 {
|
Chris@49
|
33 typedef float T;
|
Chris@49
|
34 return eT( arma_atlas(cblas_sdot)(N, (const T*)X, 1, (const T*)Y, 1) );
|
Chris@49
|
35 }
|
Chris@49
|
36 else
|
Chris@49
|
37 if(is_double<eT>::value == true)
|
Chris@49
|
38 {
|
Chris@49
|
39 typedef double T;
|
Chris@49
|
40 return eT( arma_atlas(cblas_ddot)(N, (const T*)X, 1, (const T*)Y, 1) );
|
Chris@49
|
41 }
|
Chris@49
|
42 else
|
Chris@49
|
43 {
|
Chris@49
|
44 return eT(0);
|
Chris@49
|
45 }
|
Chris@49
|
46 }
|
Chris@49
|
47
|
Chris@49
|
48
|
Chris@49
|
49
|
Chris@49
|
50 template<typename eT>
|
Chris@49
|
51 arma_inline
|
Chris@49
|
52 eT
|
Chris@49
|
53 cx_cblas_dot(const int N, const eT* X, const eT* Y)
|
Chris@49
|
54 {
|
Chris@49
|
55 arma_type_check((is_supported_blas_type<eT>::value == false));
|
Chris@49
|
56
|
Chris@49
|
57 if(is_supported_complex_float<eT>::value == true)
|
Chris@49
|
58 {
|
Chris@49
|
59 typedef typename std::complex<float> T;
|
Chris@49
|
60
|
Chris@49
|
61 T out;
|
Chris@49
|
62 arma_atlas(cblas_cdotu_sub)(N, (const T*)X, 1, (const T*)Y, 1, &out);
|
Chris@49
|
63
|
Chris@49
|
64 return eT(out);
|
Chris@49
|
65 }
|
Chris@49
|
66 else
|
Chris@49
|
67 if(is_supported_complex_double<eT>::value == true)
|
Chris@49
|
68 {
|
Chris@49
|
69 typedef typename std::complex<double> T;
|
Chris@49
|
70
|
Chris@49
|
71 T out;
|
Chris@49
|
72 arma_atlas(cblas_zdotu_sub)(N, (const T*)X, 1, (const T*)Y, 1, &out);
|
Chris@49
|
73
|
Chris@49
|
74 return eT(out);
|
Chris@49
|
75 }
|
Chris@49
|
76 else
|
Chris@49
|
77 {
|
Chris@49
|
78 return eT(0);
|
Chris@49
|
79 }
|
Chris@49
|
80 }
|
Chris@49
|
81
|
Chris@49
|
82
|
Chris@49
|
83
|
Chris@49
|
84 template<typename eT>
|
Chris@49
|
85 inline
|
Chris@49
|
86 void
|
Chris@49
|
87 cblas_gemv
|
Chris@49
|
88 (
|
Chris@49
|
89 const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA,
|
Chris@49
|
90 const int M, const int N,
|
Chris@49
|
91 const eT alpha,
|
Chris@49
|
92 const eT *A, const int lda,
|
Chris@49
|
93 const eT *X, const int incX,
|
Chris@49
|
94 const eT beta,
|
Chris@49
|
95 eT *Y, const int incY
|
Chris@49
|
96 )
|
Chris@49
|
97 {
|
Chris@49
|
98 arma_type_check((is_supported_blas_type<eT>::value == false));
|
Chris@49
|
99
|
Chris@49
|
100 if(is_float<eT>::value == true)
|
Chris@49
|
101 {
|
Chris@49
|
102 typedef float T;
|
Chris@49
|
103 arma_atlas(cblas_sgemv)(Order, TransA, M, N, (const T)tmp_real(alpha), (const T*)A, lda, (const T*)X, incX, (const T)tmp_real(beta), (T*)Y, incY);
|
Chris@49
|
104 }
|
Chris@49
|
105 else
|
Chris@49
|
106 if(is_double<eT>::value == true)
|
Chris@49
|
107 {
|
Chris@49
|
108 typedef double T;
|
Chris@49
|
109 arma_atlas(cblas_dgemv)(Order, TransA, M, N, (const T)tmp_real(alpha), (const T*)A, lda, (const T*)X, incX, (const T)tmp_real(beta), (T*)Y, incY);
|
Chris@49
|
110 }
|
Chris@49
|
111 else
|
Chris@49
|
112 if(is_supported_complex_float<eT>::value == true)
|
Chris@49
|
113 {
|
Chris@49
|
114 typedef std::complex<float> T;
|
Chris@49
|
115 arma_atlas(cblas_cgemv)(Order, TransA, M, N, (const T*)&alpha, (const T*)A, lda, (const T*)X, incX, (const T*)&beta, (T*)Y, incY);
|
Chris@49
|
116 }
|
Chris@49
|
117 else
|
Chris@49
|
118 if(is_supported_complex_double<eT>::value == true)
|
Chris@49
|
119 {
|
Chris@49
|
120 typedef std::complex<double> T;
|
Chris@49
|
121 arma_atlas(cblas_zgemv)(Order, TransA, M, N, (const T*)&alpha, (const T*)A, lda, (const T*)X, incX, (const T*)&beta, (T*)Y, incY);
|
Chris@49
|
122 }
|
Chris@49
|
123 }
|
Chris@49
|
124
|
Chris@49
|
125
|
Chris@49
|
126
|
Chris@49
|
127 template<typename eT>
|
Chris@49
|
128 inline
|
Chris@49
|
129 void
|
Chris@49
|
130 cblas_gemm
|
Chris@49
|
131 (
|
Chris@49
|
132 const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA,
|
Chris@49
|
133 const enum CBLAS_TRANSPOSE TransB, const int M, const int N,
|
Chris@49
|
134 const int K, const eT alpha, const eT *A,
|
Chris@49
|
135 const int lda, const eT *B, const int ldb,
|
Chris@49
|
136 const eT beta, eT *C, const int ldc
|
Chris@49
|
137 )
|
Chris@49
|
138 {
|
Chris@49
|
139 arma_type_check((is_supported_blas_type<eT>::value == false));
|
Chris@49
|
140
|
Chris@49
|
141 if(is_float<eT>::value == true)
|
Chris@49
|
142 {
|
Chris@49
|
143 typedef float T;
|
Chris@49
|
144 arma_atlas(cblas_sgemm)(Order, TransA, TransB, M, N, K, (const T)tmp_real(alpha), (const T*)A, lda, (const T*)B, ldb, (const T)tmp_real(beta), (T*)C, ldc);
|
Chris@49
|
145 }
|
Chris@49
|
146 else
|
Chris@49
|
147 if(is_double<eT>::value == true)
|
Chris@49
|
148 {
|
Chris@49
|
149 typedef double T;
|
Chris@49
|
150 arma_atlas(cblas_dgemm)(Order, TransA, TransB, M, N, K, (const T)tmp_real(alpha), (const T*)A, lda, (const T*)B, ldb, (const T)tmp_real(beta), (T*)C, ldc);
|
Chris@49
|
151 }
|
Chris@49
|
152 else
|
Chris@49
|
153 if(is_supported_complex_float<eT>::value == true)
|
Chris@49
|
154 {
|
Chris@49
|
155 typedef std::complex<float> T;
|
Chris@49
|
156 arma_atlas(cblas_cgemm)(Order, TransA, TransB, M, N, K, (const T*)&alpha, (const T*)A, lda, (const T*)B, ldb, (const T*)&beta, (T*)C, ldc);
|
Chris@49
|
157 }
|
Chris@49
|
158 else
|
Chris@49
|
159 if(is_supported_complex_double<eT>::value == true)
|
Chris@49
|
160 {
|
Chris@49
|
161 typedef std::complex<double> T;
|
Chris@49
|
162 arma_atlas(cblas_zgemm)(Order, TransA, TransB, M, N, K, (const T*)&alpha, (const T*)A, lda, (const T*)B, ldb, (const T*)&beta, (T*)C, ldc);
|
Chris@49
|
163 }
|
Chris@49
|
164 }
|
Chris@49
|
165
|
Chris@49
|
166
|
Chris@49
|
167
|
Chris@49
|
168 template<typename eT>
|
Chris@49
|
169 inline
|
Chris@49
|
170 int
|
Chris@49
|
171 clapack_getrf
|
Chris@49
|
172 (
|
Chris@49
|
173 const enum CBLAS_ORDER Order, const int M, const int N,
|
Chris@49
|
174 eT *A, const int lda, int *ipiv
|
Chris@49
|
175 )
|
Chris@49
|
176 {
|
Chris@49
|
177 arma_type_check((is_supported_blas_type<eT>::value == false));
|
Chris@49
|
178
|
Chris@49
|
179 if(is_float<eT>::value == true)
|
Chris@49
|
180 {
|
Chris@49
|
181 typedef float T;
|
Chris@49
|
182 return arma_atlas(clapack_sgetrf)(Order, M, N, (T*)A, lda, ipiv);
|
Chris@49
|
183 }
|
Chris@49
|
184 else
|
Chris@49
|
185 if(is_double<eT>::value == true)
|
Chris@49
|
186 {
|
Chris@49
|
187 typedef double T;
|
Chris@49
|
188 return arma_atlas(clapack_dgetrf)(Order, M, N, (T*)A, lda, ipiv);
|
Chris@49
|
189 }
|
Chris@49
|
190 else
|
Chris@49
|
191 if(is_supported_complex_float<eT>::value == true)
|
Chris@49
|
192 {
|
Chris@49
|
193 typedef std::complex<float> T;
|
Chris@49
|
194 return arma_atlas(clapack_cgetrf)(Order, M, N, (T*)A, lda, ipiv);
|
Chris@49
|
195 }
|
Chris@49
|
196 else
|
Chris@49
|
197 if(is_supported_complex_double<eT>::value == true)
|
Chris@49
|
198 {
|
Chris@49
|
199 typedef std::complex<double> T;
|
Chris@49
|
200 return arma_atlas(clapack_zgetrf)(Order, M, N, (T*)A, lda, ipiv);
|
Chris@49
|
201 }
|
Chris@49
|
202 else
|
Chris@49
|
203 {
|
Chris@49
|
204 return -1;
|
Chris@49
|
205 }
|
Chris@49
|
206 }
|
Chris@49
|
207
|
Chris@49
|
208
|
Chris@49
|
209
|
Chris@49
|
210 template<typename eT>
|
Chris@49
|
211 inline
|
Chris@49
|
212 int
|
Chris@49
|
213 clapack_getri
|
Chris@49
|
214 (
|
Chris@49
|
215 const enum CBLAS_ORDER Order, const int N, eT *A,
|
Chris@49
|
216 const int lda, const int *ipiv
|
Chris@49
|
217 )
|
Chris@49
|
218 {
|
Chris@49
|
219 arma_type_check((is_supported_blas_type<eT>::value == false));
|
Chris@49
|
220
|
Chris@49
|
221 if(is_float<eT>::value == true)
|
Chris@49
|
222 {
|
Chris@49
|
223 typedef float T;
|
Chris@49
|
224 return arma_atlas(clapack_sgetri)(Order, N, (T*)A, lda, ipiv);
|
Chris@49
|
225 }
|
Chris@49
|
226 else
|
Chris@49
|
227 if(is_double<eT>::value == true)
|
Chris@49
|
228 {
|
Chris@49
|
229 typedef double T;
|
Chris@49
|
230 return arma_atlas(clapack_dgetri)(Order, N, (T*)A, lda, ipiv);
|
Chris@49
|
231 }
|
Chris@49
|
232 else
|
Chris@49
|
233 if(is_supported_complex_float<eT>::value == true)
|
Chris@49
|
234 {
|
Chris@49
|
235 typedef std::complex<float> T;
|
Chris@49
|
236 return arma_atlas(clapack_cgetri)(Order, N, (T*)A, lda, ipiv);
|
Chris@49
|
237 }
|
Chris@49
|
238 else
|
Chris@49
|
239 if(is_supported_complex_double<eT>::value == true)
|
Chris@49
|
240 {
|
Chris@49
|
241 typedef std::complex<double> T;
|
Chris@49
|
242 return arma_atlas(clapack_zgetri)(Order, N, (T*)A, lda, ipiv);
|
Chris@49
|
243 }
|
Chris@49
|
244 else
|
Chris@49
|
245 {
|
Chris@49
|
246 return -1;
|
Chris@49
|
247 }
|
Chris@49
|
248 }
|
Chris@49
|
249
|
Chris@49
|
250
|
Chris@49
|
251
|
Chris@49
|
252 template<typename eT>
|
Chris@49
|
253 inline
|
Chris@49
|
254 int
|
Chris@49
|
255 clapack_gesv
|
Chris@49
|
256 (
|
Chris@49
|
257 const enum CBLAS_ORDER Order,
|
Chris@49
|
258 const int N, const int NRHS,
|
Chris@49
|
259 eT* A, const int lda, int* ipiv,
|
Chris@49
|
260 eT* B, const int ldb
|
Chris@49
|
261 )
|
Chris@49
|
262 {
|
Chris@49
|
263 arma_type_check((is_supported_blas_type<eT>::value == false));
|
Chris@49
|
264
|
Chris@49
|
265 if(is_float<eT>::value == true)
|
Chris@49
|
266 {
|
Chris@49
|
267 typedef float T;
|
Chris@49
|
268 return arma_atlas(clapack_sgesv)(Order, N, NRHS, (T*)A, lda, ipiv, (T*)B, ldb);
|
Chris@49
|
269 }
|
Chris@49
|
270 else
|
Chris@49
|
271 if(is_double<eT>::value == true)
|
Chris@49
|
272 {
|
Chris@49
|
273 typedef double T;
|
Chris@49
|
274 return arma_atlas(clapack_dgesv)(Order, N, NRHS, (T*)A, lda, ipiv, (T*)B, ldb);
|
Chris@49
|
275 }
|
Chris@49
|
276 else
|
Chris@49
|
277 if(is_supported_complex_float<eT>::value == true)
|
Chris@49
|
278 {
|
Chris@49
|
279 typedef std::complex<float> T;
|
Chris@49
|
280 return arma_atlas(clapack_cgesv)(Order, N, NRHS, (T*)A, lda, ipiv, (T*)B, ldb);
|
Chris@49
|
281 }
|
Chris@49
|
282 else
|
Chris@49
|
283 if(is_supported_complex_double<eT>::value == true)
|
Chris@49
|
284 {
|
Chris@49
|
285 typedef std::complex<double> T;
|
Chris@49
|
286 return arma_atlas(clapack_zgesv)(Order, N, NRHS, (T*)A, lda, ipiv, (T*)B, ldb);
|
Chris@49
|
287 }
|
Chris@49
|
288 else
|
Chris@49
|
289 {
|
Chris@49
|
290 return -1;
|
Chris@49
|
291 }
|
Chris@49
|
292 }
|
Chris@49
|
293
|
Chris@49
|
294
|
Chris@49
|
295
|
Chris@49
|
296 }
|
Chris@49
|
297
|
Chris@49
|
298 #endif
|