[boost] Re: Matrix Multiplication Benchmark: uBLAS, Matlab 6.5, Matlab 5.3, ATLAS

From: Patrick Kowalzick (Patrick.Kowalzick_at_[hidden])
Date: 2003-05-16 02:47:13

Hello all,

I just want to send you a small table with time needed for multiplications
in several cases. The first column represent the used data_types
R=row_major,C=column_major. The 6 last columns are classical multiplications
with loops. rkc->rck is just a differnent order of the loops. I used uBlas
there as the container and I think the result is pretty close to the use of
valarray or double*.

Big surprise: Take the right configuration and uBlas runs extremly!!! fast,
but taking the wrong one.....
I will try it for 1500x1500 now as well, but as far as I have seen, it will
produce more or less the same result.

Patrick

size of metrices - 500x500
      ublas     rkc     rck     krc     kcr     crk     ckr
RRR   2.453   2.218   4.141   2.515   5.344   3.343   4.937
RRC   1.093    4.36   2.329   3.235   5.375   2.359   4.844
RCR   3.922   2.391   4.984   2.437     4.5   4.782   2.859
RCC    1.86   4.235   3.406    3.14   4.437    4.25   2.813
CRR   1.906   2.875   4.203   4.422   3.282    3.25    4.25
CRC   1.078   5.109   2.453   5.437   3.125   2.407   4.172
CCR   3.922   2.937   4.875   4.531   2.563   4.828   2.422
CCC   2.438    4.89   2.984    5.25   2.516   4.062   2.375

YES, I FEEL SORRY TO POST CODE SNIPPETS WITH THESE MACROS IN THE C++
COMMUNITY, ... anyway.

**********************code used*************************
#include <iostream>
#include <iomanip>

#include <boost/numeric/ublas/matrix.hpp>
#include <boost/numeric/ublas/io.hpp>
#include <boost/timer.hpp>

namespace ublas = boost::numeric::ublas;
using namespace std;

void initmatrix(ublas::matrix<double>& A) {
 for(size_t i = 0; i<A.size1(); i++)
  for(size_t j = 0; j<A.size2(); j++)
   A(i,j) = 1.0*rand();
}
void initmatrix(ublas::matrix<double,ublas::column_major>& A) {
 for(size_t i = 0; i<A.size1(); i++)
  for(size_t j = 0; j<A.size2(); j++)
   A(i,j) = 1.0*rand();
}

typedef ublas::matrix<double,ublas::column_major> CT;
typedef ublas::matrix<double,ublas::row_major> RT;

#define mul(a,b,d) \
{ \
  initmatrix(A); \
  initmatrix(B); \
  t.restart(); \
  for (size_t a=0;a<size;a++) \
   for (size_t b=0;b<size;b++) \
    for (size_t d=0;d<size;d++) \
     X(r,c) += A(r,k)*B(k,c); \
  time = t.elapsed(); \
  cout << setw(8) << time; \
}

#define mul2(xtype) \
{ \
  initmatrix(A); \
  initmatrix(B); \
  t.restart(); \
  X = ublas::prod< xtype >(A,B); \
  time = t.elapsed(); \
  cout << setw(8) << time; \
}

#define mulset(xtype,atype,btype) \
{ \
 xtype X(size,size); \
 atype A(size,size); \
 btype B(size,size); \
 mul2(xtype); \
 mul(r,k,c); \
 mul(r,c,k); \
 mul(k,r,c); \
 mul(k,c,r); \
 mul(c,r,k); \
 mul(c,k,r); \
 cout << endl; \
}

void test2(const size_t size) {
  boost::timer t;
        double time;
  cout << "size of metrices - " << size << "x" << size << endl << endl;
  cout << "      ublas     rkc     rck     krc     kcr     crk     ckr" <<
endl;
  cout << "RRR";mulset(RT,RT,RT);
  cout << "RRC";mulset(RT,RT,CT);
  cout << "RCR";mulset(RT,CT,RT);
  cout << "RCC";mulset(RT,CT,CT);
  //cout << "CRR";mulset(CT,RT,RT);
  //cout << "CRC";mulset(CT,RT,CT);
  //cout << "CCR";mulset(CT,CT,RT);
  //cout << "CCC";mulset(CT,CT,CT);
}

void main() {
  test2(500);
}