32 # include <itpp/config.h>
34 # include <itpp/config_msvc.h>
37 #if defined (HAVE_BLAS)
38 # include <itpp/base/blas.h>
49 cmat temp(no_cols, no_rows);
50 for (
int i = 0; i < no_rows; i++)
51 for (
int j = 0; j < no_cols; j++)
60 #if defined(HAVE_BLAS)
63 mat& mat::operator*=(
const mat &m)
65 it_assert_debug(no_cols == m.no_rows,
"mat::operator*=(): Wrong sizes");
66 mat r(no_rows, m.no_cols);
70 blas::dgemm_(&trans, &trans, &no_rows, &m.no_cols, &no_cols, &alpha, data,
71 &no_rows, m.data, &m.no_rows, &beta, r.data, &r.no_rows);
77 cmat& cmat::operator*=(
const cmat &m)
79 it_assert_debug(no_cols == m.no_rows,
"cmat::operator*=(): Wrong sizes");
80 cmat r(no_rows, m.no_cols);
81 std::complex<double> alpha = std::complex<double>(1.0);
82 std::complex<double> beta = std::complex<double>(0.0);
84 blas::zgemm_(&trans, &trans, &no_rows, &m.no_cols, &no_cols, &alpha, data,
85 &no_rows, m.data, &m.no_rows, &beta, r.data, &r.no_rows);
91 mat& mat::operator*=(
const mat &m)
93 it_assert_debug(no_cols == m.no_rows,
"Mat<>::operator*=(): Wrong sizes");
94 mat r(no_rows, m.no_cols);
95 int r_pos = 0, pos = 0, m_pos = 0;
97 for (
int i = 0; i < r.no_cols; i++) {
98 for (
int j = 0; j < r.no_rows; j++) {
101 for (
int k = 0; k < no_cols; k++) {
102 tmp += data[pos+j] * m.data[m_pos+k];
105 r.data[r_pos+j] = tmp;
115 cmat& cmat::operator*=(
const cmat &m)
117 it_assert_debug(no_cols == m.no_rows,
"Mat<>::operator*=(): Wrong sizes");
118 cmat r(no_rows, m.no_cols);
119 int r_pos = 0, pos = 0, m_pos = 0;
121 for (
int i = 0; i < r.no_cols; i++) {
122 for (
int j = 0; j < r.no_rows; j++) {
123 std::complex<double> tmp(0.0);
125 for (
int k = 0; k < no_cols; k++) {
126 tmp += data[pos+j] * m.data[m_pos+k];
129 r.data[r_pos+j] = tmp;
140 #if defined(HAVE_BLAS)
142 mat
operator*(
const mat &m1,
const mat &m2)
144 it_assert_debug(m1.no_cols == m2.no_rows,
"mat::operator*(): Wrong sizes");
145 mat r(m1.no_rows, m2.no_cols);
149 blas::dgemm_(&trans, &trans, &m1.no_rows, &m2.no_cols, &m1.no_cols, &alpha,
150 m1.data, &m1.no_rows, m2.data, &m2.no_rows, &beta, r.data,
156 cmat
operator*(
const cmat &m1,
const cmat &m2)
158 it_assert_debug(m1.no_cols == m2.no_rows,
"cmat::operator*(): Wrong sizes");
159 cmat r(m1.no_rows, m2.no_cols);
160 std::complex<double> alpha = std::complex<double>(1.0);
161 std::complex<double> beta = std::complex<double>(0.0);
163 blas::zgemm_(&trans, &trans, &m1.no_rows, &m2.no_cols, &m1.no_cols, &alpha,
164 m1.data, &m1.no_rows, m2.data, &m2.no_rows, &beta, r.data,
170 mat
operator*(
const mat &m1,
const mat &m2)
173 "Mat<>::operator*(): Wrong sizes");
174 mat r(m1.no_rows, m2.no_cols);
177 double *t2 = m2.data;
178 for (
int i = 0; i < r.no_cols; i++) {
179 for (
int j = 0; j < r.no_rows; j++) {
182 for (
int k = m1.no_cols; k > 0; k--) {
183 tmp += *(t1) * *(t2++);
195 cmat
operator*(
const cmat &m1,
const cmat &m2)
198 "Mat<>::operator*(): Wrong sizes");
199 cmat r(m1.no_rows, m2.no_cols);
200 std::complex<double> *tr = r.data;
201 std::complex<double> *t1;
202 std::complex<double> *t2 = m2.data;
203 for (
int i = 0; i < r.no_cols; i++) {
204 for (
int j = 0; j < r.no_rows; j++) {
205 std::complex<double> tmp(0.0);
207 for (
int k = m1.no_cols; k > 0; k--) {
208 tmp += *(t1) * *(t2++);
221 #if defined(HAVE_BLAS)
223 vec
operator*(
const mat &m,
const vec &v)
225 it_assert_debug(m.no_cols == v.size(),
"mat::operator*(): Wrong sizes");
231 blas::dgemv_(&trans, &m.no_rows, &m.no_cols, &alpha, m.data, &m.no_rows,
232 v._data(), &incr, &beta, r._data(), &incr);
237 cvec
operator*(
const cmat &m,
const cvec &v)
239 it_assert_debug(m.no_cols == v.size(),
"cmat::operator*(): Wrong sizes");
241 std::complex<double> alpha = std::complex<double>(1.0);
242 std::complex<double> beta = std::complex<double>(0.0);
245 blas::zgemv_(&trans, &m.no_rows, &m.no_cols, &alpha, m.data, &m.no_rows,
246 v._data(), &incr, &beta, r._data(), &incr);
251 vec
operator*(
const mat &m,
const vec &v)
254 "Mat<>::operator*(): Wrong sizes");
256 for (
int i = 0; i < m.no_rows; i++) {
259 for (
int k = 0; k < m.no_cols; k++) {
260 r(i) += m.data[m_pos+i] * v(k);
268 cvec
operator*(
const cmat &m,
const cvec &v)
271 "Mat<>::operator*(): Wrong sizes");
273 for (
int i = 0; i < m.no_rows; i++) {
274 r(i) = std::complex<double>(0.0);
276 for (
int k = 0; k < m.no_cols; k++) {
277 r(i) += m.data[m_pos+i] * v(k);
292 template class Mat<double>;
293 template class Mat<std::complex<double> >;
294 template class Mat<int>;
295 template class Mat<short int>;
296 template class Mat<bin>;
300 template mat
operator+(
const mat &m1,
const mat &m2);
301 template cmat
operator+(
const cmat &m1,
const cmat &m2);
302 template imat
operator+(
const imat &m1,
const imat &m2);
303 template smat
operator+(
const smat &m1,
const smat &m2);
306 template mat
operator+(
const mat &m,
double t);
307 template cmat
operator+(
const cmat &m, std::complex<double> t);
308 template imat
operator+(
const imat &m,
int t);
309 template smat
operator+(
const smat &m,
short t);
312 template mat
operator+(
double t,
const mat &m);
313 template cmat
operator+(std::complex<double> t,
const cmat &m);
314 template imat
operator+(
int t,
const imat &m);
315 template smat
operator+(
short t,
const smat &m);
320 template mat
operator-(
const mat &m1,
const mat &m2);
321 template cmat
operator-(
const cmat &m1,
const cmat &m2);
322 template imat
operator-(
const imat &m1,
const imat &m2);
323 template smat
operator-(
const smat &m1,
const smat &m2);
326 template mat
operator-(
const mat &m,
double t);
327 template cmat
operator-(
const cmat &m, std::complex<double> t);
328 template imat
operator-(
const imat &m,
int t);
329 template smat
operator-(
const smat &m,
short t);
332 template mat
operator-(
double t,
const mat &m);
333 template cmat
operator-(std::complex<double> t,
const cmat &m);
334 template imat
operator-(
int t,
const imat &m);
335 template smat
operator-(
short t,
const smat &m);
348 template imat
operator*(
const imat &m1,
const imat &m2);
349 template smat
operator*(
const smat &m1,
const smat &m2);
352 template ivec
operator*(
const imat &m,
const ivec &v);
353 template svec
operator*(
const smat &m,
const svec &v);
356 template mat
operator*(
const mat &m,
double t);
357 template cmat
operator*(
const cmat &m, std::complex<double> t);
358 template imat
operator*(
const imat &m,
int t);
359 template smat
operator*(
const smat &m,
short t);
362 template mat
operator*(
double t,
const mat &m);
363 template cmat
operator*(std::complex<double> t,
const cmat &m);
364 template imat
operator*(
int t,
const imat &m);
365 template smat
operator*(
short t,
const smat &m);
370 template mat
elem_mult(
const mat &m1,
const mat &m2);
371 template cmat
elem_mult(
const cmat &m1,
const cmat &m2);
372 template imat
elem_mult(
const imat &m1,
const imat &m2);
373 template smat
elem_mult(
const smat &m1,
const smat &m2);
376 template void elem_mult_out(
const mat &m1,
const mat &m2, mat &out);
377 template void elem_mult_out(
const cmat &m1,
const cmat &m2, cmat &out);
378 template void elem_mult_out(
const imat &m1,
const imat &m2, imat &out);
379 template void elem_mult_out(
const smat &m1,
const smat &m2, smat &out);
383 const mat &m3, mat &out);
385 const cmat &m3, cmat &out);
387 const imat &m3, imat &out);
389 const smat &m3, smat &out);
393 template void elem_mult_out(
const mat &m1,
const mat &m2,
const mat &m3,
394 const mat &m4, mat &out);
396 const cmat &m3,
const cmat &m4, cmat &out);
398 const imat &m3,
const imat &m4, imat &out);
400 const smat &m3,
const smat &m4, smat &out);
411 template std::complex<double>
elem_mult_sum(
const cmat &m1,
const cmat &m2);
413 template short elem_mult_sum(
const smat &m1,
const smat &m2);
418 template mat
operator/(
double t,
const mat &m);
419 template cmat
operator/(std::complex<double> t,
const cmat &m);
420 template imat
operator/(
int t,
const imat &m);
421 template smat
operator/(
short t,
const smat &m);
424 template mat
operator/(
const mat &m,
double t);
425 template cmat
operator/(
const cmat &m, std::complex<double> t);
426 template imat
operator/(
const imat &m,
int t);
427 template smat
operator/(
const smat &m,
short t);
432 template mat
elem_div(
const mat &m1,
const mat &m2);
433 template cmat
elem_div(
const cmat &m1,
const cmat &m2);
434 template imat
elem_div(
const imat &m1,
const imat &m2);
435 template smat
elem_div(
const smat &m1,
const smat &m2);
438 template void elem_div_out(
const mat &m1,
const mat &m2, mat &out);
439 template void elem_div_out(
const cmat &m1,
const cmat &m2, cmat &out);
440 template void elem_div_out(
const imat &m1,
const imat &m2, imat &out);
441 template void elem_div_out(
const smat &m1,
const smat &m2, smat &out);
444 template double elem_div_sum(
const mat &m1,
const mat &m2);
445 template std::complex<double>
elem_div_sum(
const cmat &m1,
447 template int elem_div_sum(
const imat &m1,
const imat &m2);
448 template short elem_div_sum(
const smat &m1,
const smat &m2);
467 template std::ostream &
operator<<(std::ostream &os,
const mat &m);
468 template std::ostream &
operator<<(std::ostream &os,
const cmat &m);
469 template std::ostream &
operator<<(std::ostream &os,
const imat &m);
470 template std::ostream &
operator<<(std::ostream &os,
const smat &m);
471 template std::ostream &
operator<<(std::ostream &os,
const bmat &m);
473 template std::istream &
operator>>(std::istream &is, mat &m);
474 template std::istream &
operator>>(std::istream &is, cmat &m);
475 template std::istream &
operator>>(std::istream &is, imat &m);
476 template std::istream &
operator>>(std::istream &is, smat &m);