#include <blaze/Math.h> #include <iostream> #include <random> #include "gemm.h" // fill rectangular matrix with random values template <typename MATRIX> void fill(MATRIX &A) { typedef typename MATRIX::ElementType T; std::random_device random; std::default_random_engine mt(random()); std::uniform_real_distribution<T> uniform(-100,100); for (std::size_t i=0; i<(~A).rows(); ++i) { for (std::size_t j=0; j<(~A).columns(); ++j) { A(i,j) = uniform(mt); } } } int main() { typedef double ElementType; constexpr auto SOA = blaze::rowMajor; constexpr auto SOB = blaze::rowMajor; constexpr auto SOC = blaze::rowMajor; std::size_t m = 5; std::size_t n = 5; std::size_t k = 5; ElementType alpha = 1; ElementType beta = 0; blaze::DynamicMatrix<double, SOA> A(m, k); blaze::DynamicMatrix<double, SOB> B(k, n); blaze::DynamicMatrix<double, SOC> C(m, n); fill(A); fill(B); fill(C); std::cout << "A = " << A << std::endl; std::cout << "B = " << B << std::endl; std::cout << "C = " << C << std::endl; std::cout << "C <- " << beta << "*C + " << alpha << "*A*B" << std::endl; foo::gemm(alpha, A, B, beta, C); std::cout << "C = " << C << std::endl; std::cout << "C <- " << beta << "*C + " << alpha << "*(A^T+A)*(2*B)" << std::endl; foo::gemm(alpha, blaze::trans(A)+A, 2*B, beta, C); std::cout << "C = " << C << std::endl; } |