| 1 |  | 
|---|
| 2 | #include <itpp/itbase.h> | 
|---|
| 3 | using namespace itpp; | 
|---|
| 4 |  | 
|---|
| 5 | //These lines are needed for use of cout and endl | 
|---|
| 6 | using std::cout; | 
|---|
| 7 | using std::endl; | 
|---|
| 8 |  | 
|---|
| 9 | //#define NDEBUG | 
|---|
| 10 |  | 
|---|
| 11 | mat matmul ( mat &A, mat &B ) { | 
|---|
| 12 |         mat C ( A.rows(),B.cols() ); | 
|---|
| 13 |         int i,j,k; | 
|---|
| 14 |         double sum; | 
|---|
| 15 |  | 
|---|
| 16 |         for ( i=0;i<A.rows();i++ ) { | 
|---|
| 17 |                 for ( j=0;j<A.cols();j++ ) { | 
|---|
| 18 |                         sum = 0.0; | 
|---|
| 19 |                         for ( k=0;k<A.cols();k++ ) { | 
|---|
| 20 |                                 sum+=A._elem ( i,k ) *B._elem ( k,j ); | 
|---|
| 21 |                         } | 
|---|
| 22 |                         C ( i,j ) = sum; | 
|---|
| 23 |                 } | 
|---|
| 24 |         } | 
|---|
| 25 |         return C; | 
|---|
| 26 | } | 
|---|
| 27 |  | 
|---|
| 28 | void matmul2 ( int n,  double *A, double *B, double *C ) { | 
|---|
| 29 |         int i,j,k; | 
|---|
| 30 |         double sum; | 
|---|
| 31 |  | 
|---|
| 32 |         for ( i=0;i<n;i++ ) { | 
|---|
| 33 |                 for ( j=0;j<n;j++ ) { | 
|---|
| 34 |                         sum = 0.0; | 
|---|
| 35 |                         for ( k=0;k<n;k++ ) { | 
|---|
| 36 |                                 sum+=A [ i*n+k ] * B [ k*n+j ]; | 
|---|
| 37 |                         } | 
|---|
| 38 |                         C[ i*n+j] = sum; | 
|---|
| 39 |                 } | 
|---|
| 40 |         } | 
|---|
| 41 | //      return C; | 
|---|
| 42 | } | 
|---|
| 43 |  | 
|---|
| 44 | int main() { | 
|---|
| 45 |         Real_Timer tt; | 
|---|
| 46 |         vec exec_times ( 4 ); | 
|---|
| 47 |         vec exec_times_b ( 4 ); | 
|---|
| 48 |         vec exec_times_c ( 4 ); | 
|---|
| 49 |  | 
|---|
| 50 |         mat A; | 
|---|
| 51 |         mat B; | 
|---|
| 52 |         mat C; | 
|---|
| 53 |  | 
|---|
| 54 |         vec vn="5 50 200 500"; | 
|---|
| 55 |         int n; | 
|---|
| 56 |  | 
|---|
| 57 |         for ( int i=0;i<vn.length();i++ ) { | 
|---|
| 58 |                 n = vn ( i ); | 
|---|
| 59 |                 A = randu ( n,n ); | 
|---|
| 60 |                 B = randu ( n,n ); | 
|---|
| 61 |  | 
|---|
| 62 |                 tt.tic(); | 
|---|
| 63 |                 for ( int ii=0;ii<10;ii++ ) {C = matmul ( A,B );} | 
|---|
| 64 |                 exec_times ( i ) =tt.toc(); | 
|---|
| 65 |  | 
|---|
| 66 |                 tt.tic(); | 
|---|
| 67 |                 for ( int ii=0;ii<10;ii++ ) {C = A*B;} | 
|---|
| 68 |                 exec_times_b ( i ) =tt.toc(); | 
|---|
| 69 |  | 
|---|
| 70 |                 C = zeros(n,n); | 
|---|
| 71 |                 tt.tic(); | 
|---|
| 72 |                 for ( int ii=0;ii<10;ii++ ) { matmul2(n,A._data(),B._data(),C._data());} | 
|---|
| 73 |                 exec_times_c ( i ) =tt.toc(); | 
|---|
| 74 |         } | 
|---|
| 75 |         cout << exec_times <<endl; | 
|---|
| 76 |         cout << exec_times_b <<endl; | 
|---|
| 77 |         cout << exec_times_c <<endl; | 
|---|
| 78 |  | 
|---|
| 79 |         it_file itf ( "blas_test.it" ); | 
|---|
| 80 |         itf << Name ( "exec_times" ) <<exec_times; | 
|---|
| 81 |         itf << Name ( "exec_times_b" ) <<exec_times_b; | 
|---|
| 82 |         itf << Name ( "exec_times_c" ) <<exec_times_c; | 
|---|
| 83 |  | 
|---|
| 84 |         return 0; | 
|---|
| 85 | } | 
|---|