1 //
2 //  Copyright (c) 2000-2002
3 //  Joerg Walter, Mathias Koch
4 //
5 //  Distributed under the Boost Software License, Version 1.0. (See
6 //  accompanying file LICENSE_1_0.txt or copy at
7 //  http://www.boost.org/LICENSE_1_0.txt)
8 //
9 //  The authors gratefully acknowledge the support of
10 //  GeNeSys mbH & Co. KG in producing this work.
11 //
12 
13 #include "bench2.hpp"
14 
15 template<class T, int N>
16 struct bench_c_outer_prod {
17     typedef T value_type;
18 
operator ()bench_c_outer_prod19     void operator () (int runs) const {
20         try {
21             static typename c_matrix_traits<T, N, N>::type m;
22             static typename c_vector_traits<T, N>::type v1, v2;
23             initialize_c_vector<T, N> () (v1);
24             initialize_c_vector<T, N> () (v2);
25             boost::timer t;
26             for (int i = 0; i < runs; ++ i) {
27                 for (int j = 0; j < N; ++ j) {
28                     for (int k = 0; k < N; ++ k) {
29                         m [j] [k] = - v1 [j] * v2 [k];
30                     }
31                 }
32 //                sink_c_matrix<T, N, N> () (m);
33                 BOOST_UBLAS_NOT_USED(m);
34             }
35             footer<value_type> () (N * N, N * N, runs, t.elapsed ());
36         }
37         catch (std::exception &e) {
38             std::cout << e.what () << std::endl;
39         }
40     }
41 };
42 template<class M, class V, int N>
43 struct bench_my_outer_prod {
44     typedef typename M::value_type value_type;
45 
operator ()bench_my_outer_prod46     void operator () (int runs, safe_tag) const {
47         try {
48             static M m (N, N, N * N);
49             static V v1 (N, N), v2 (N, N);
50             initialize_vector (v1);
51             initialize_vector (v2);
52             boost::timer t;
53             for (int i = 0; i < runs; ++ i) {
54                 m = - ublas::outer_prod (v1, v2);
55 //                sink_matrix (m);
56                 BOOST_UBLAS_NOT_USED(m);
57             }
58             footer<value_type> () (N * N, N * N, runs, t.elapsed ());
59         }
60         catch (std::exception &e) {
61             std::cout << e.what () << std::endl;
62         }
63     }
operator ()bench_my_outer_prod64     void operator () (int runs, fast_tag) const {
65         try {
66             static M m (N, N, N * N);
67             static V v1 (N, N), v2 (N, N);
68             initialize_vector (v1);
69             initialize_vector (v2);
70             boost::timer t;
71             for (int i = 0; i < runs; ++ i) {
72                 m.assign (- ublas::outer_prod (v1, v2));
73 //                sink_matrix (m);
74                 BOOST_UBLAS_NOT_USED(m);
75             }
76             footer<value_type> () (N * N, N * N, runs, t.elapsed ());
77         }
78         catch (std::exception &e) {
79             std::cout << e.what () << std::endl;
80         }
81     }
82 };
83 template<class M, class V, int N>
84 struct bench_cpp_outer_prod {
85     typedef typename M::value_type value_type;
86 
operator ()bench_cpp_outer_prod87     void operator () (int runs) const {
88         try {
89             static M m (N * N);
90             static V v1 (N), v2 (N);
91             initialize_vector (v1);
92             initialize_vector (v2);
93             boost::timer t;
94             for (int i = 0; i < runs; ++ i) {
95                 for (int j = 0; j < N; ++ j) {
96                     for (int k = 0; k < N; ++ k) {
97                         m [N * j + k] = - v1 [j] * v2 [k];
98                     }
99                 }
100 //                sink_vector (m);
101             }
102             footer<value_type> () (N * N, N * N, runs, t.elapsed ());
103         }
104         catch (std::exception &e) {
105             std::cout << e.what () << std::endl;
106         }
107     }
108 };
109 
110 template<class T, int N>
111 struct bench_c_matrix_vector_prod {
112     typedef T value_type;
113 
operator ()bench_c_matrix_vector_prod114     void operator () (int runs) const {
115         try {
116             static typename c_matrix_traits<T, N, N>::type m;
117             static typename c_vector_traits<T, N>::type v1, v2;
118             initialize_c_matrix<T, N, N> () (m);
119             initialize_c_vector<T, N> () (v1);
120             boost::timer t;
121             for (int i = 0; i < runs; ++ i) {
122                 for (int j = 0; j < N; ++ j) {
123                     v2 [j] = 0;
124                     for (int k = 0; k < N; ++ k) {
125                         v2 [j] += m [j] [k] * v1 [k];
126                     }
127                 }
128 //                sink_c_vector<T, N> () (v2);
129             }
130             footer<value_type> () (N * N, N * (N - 1), runs, t.elapsed ());
131         }
132         catch (std::exception &e) {
133             std::cout << e.what () << std::endl;
134         }
135     }
136 };
137 template<class M, class V, int N>
138 struct bench_my_matrix_vector_prod {
139     typedef typename M::value_type value_type;
140 
operator ()bench_my_matrix_vector_prod141     void operator () (int runs, safe_tag) const {
142         try {
143             static M m (N, N, N * N);
144             static V v1 (N, N), v2 (N, N);
145             initialize_matrix (m);
146             initialize_vector (v1);
147             boost::timer t;
148             for (int i = 0; i < runs; ++ i) {
149                 v2 = ublas::prod (m, v1);
150 //                sink_vector (v2);
151             }
152             footer<value_type> () (N * N, N * (N - 1), runs, t.elapsed ());
153         }
154         catch (std::exception &e) {
155             std::cout << e.what () << std::endl;
156         }
157     }
operator ()bench_my_matrix_vector_prod158     void operator () (int runs, fast_tag) const {
159         try {
160             static M m (N, N, N * N);
161             static V v1 (N, N), v2 (N, N);
162             initialize_matrix (m);
163             initialize_vector (v1);
164             boost::timer t;
165             for (int i = 0; i < runs; ++ i) {
166                 v2.assign (ublas::prod (m, v1));
167 //                sink_vector (v2);
168             }
169             footer<value_type> () (N * N, N * (N - 1), runs, t.elapsed ());
170         }
171         catch (std::exception &e) {
172             std::cout << e.what () << std::endl;
173         }
174     }
175 };
176 template<class M, class V, int N>
177 struct bench_cpp_matrix_vector_prod {
178     typedef typename M::value_type value_type;
179 
operator ()bench_cpp_matrix_vector_prod180     void operator () (int runs) const {
181         try {
182             static M m (N * N);
183             static V v1 (N), v2 (N);
184             initialize_vector (m);
185             initialize_vector (v1);
186             boost::timer t;
187             for (int i = 0; i < runs; ++ i) {
188                 for (int j = 0; j < N; ++ j) {
189                     std::valarray<value_type> row (m [std::slice (N * j, N, 1)]);
190                     v2 [j] = (row * v1).sum ();
191                 }
192 //                sink_vector (v2);
193             }
194             footer<value_type> () (N * N, N * (N - 1), runs, t.elapsed ());
195         }
196         catch (std::exception &e) {
197             std::cout << e.what () << std::endl;
198         }
199     }
200 };
201 
202 template<class T, int N>
203 struct bench_c_matrix_add {
204     typedef T value_type;
205 
operator ()bench_c_matrix_add206     void operator () (int runs) const {
207         try {
208             static typename c_matrix_traits<T, N, N>::type m1, m2, m3;
209             initialize_c_matrix<T, N, N> () (m1);
210             initialize_c_matrix<T, N, N> () (m2);
211             boost::timer t;
212             for (int i = 0; i < runs; ++ i) {
213                 for (int j = 0; j < N; ++ j) {
214                     for (int k = 0; k < N; ++ k) {
215                         m3 [j] [k] = - (m1 [j] [k] + m2 [j] [k]);
216                     }
217                 }
218 //                sink_c_matrix<T, N, N> () (m3);
219                 BOOST_UBLAS_NOT_USED(m3);
220             }
221             footer<value_type> () (0, 2 * N * N, runs, t.elapsed ());
222         }
223         catch (std::exception &e) {
224             std::cout << e.what () << std::endl;
225         }
226     }
227 };
228 template<class M, int N>
229 struct bench_my_matrix_add {
230     typedef typename M::value_type value_type;
231 
operator ()bench_my_matrix_add232     void operator () (int runs, safe_tag) const {
233         try {
234             static M m1 (N, N, N * N), m2 (N, N, N * N), m3 (N, N, N * N);
235             initialize_matrix (m1);
236             initialize_matrix (m2);
237             boost::timer t;
238             for (int i = 0; i < runs; ++ i) {
239                 m3 = - (m1 + m2);
240 //                sink_matrix (m3);
241             }
242             footer<value_type> () (0, 2 * N * N, runs, t.elapsed ());
243         }
244         catch (std::exception &e) {
245             std::cout << e.what () << std::endl;
246         }
247     }
operator ()bench_my_matrix_add248     void operator () (int runs, fast_tag) const {
249         try {
250             static M m1 (N, N, N * N), m2 (N, N, N * N), m3 (N, N, N * N);
251             initialize_matrix (m1);
252             initialize_matrix (m2);
253             boost::timer t;
254             for (int i = 0; i < runs; ++ i) {
255                 m3.assign (- (m1 + m2));
256 //                sink_matrix (m3);
257             }
258             footer<value_type> () (0, 2 * N * N, runs, t.elapsed ());
259         }
260         catch (std::exception &e) {
261             std::cout << e.what () << std::endl;
262         }
263     }
264 };
265 template<class M, int N>
266 struct bench_cpp_matrix_add {
267     typedef typename M::value_type value_type;
268 
operator ()bench_cpp_matrix_add269     void operator () (int runs) const {
270         try {
271             static M m1 (N * N), m2 (N * N), m3 (N * N);
272             initialize_vector (m1);
273             initialize_vector (m2);
274             boost::timer t;
275             for (int i = 0; i < runs; ++ i) {
276                 m3 = - (m1 + m2);
277 //                sink_vector (m3);
278             }
279             footer<value_type> () (0, 2 * N * N, runs, t.elapsed ());
280         }
281         catch (std::exception &e) {
282             std::cout << e.what () << std::endl;
283         }
284     }
285 };
286 
287 // Benchmark O (n ^ 2)
288 template<class T, int N>
operator ()(int runs)289 void bench_2<T, N>::operator () (int runs) {
290     header ("bench_2");
291 
292     header ("outer_prod");
293 
294     header ("C array");
295     bench_c_outer_prod<T, N> () (runs);
296 
297 #ifdef USE_SPARSE_MATRIX
298 #ifdef USE_MAP_ARRAY
299     header ("sparse_matrix<map_array>, sparse_vector<map_array> safe");
300     bench_my_outer_prod<ublas::sparse_matrix<T, ublas::row_major, ublas::map_array<std::size_t, T> >,
301                         ublas::sparse_vector<T, ublas::map_array<std::size_t, T> >, N> () (runs, safe_tag ());
302 
303     header ("sparse_matrix<map_array>, sparse_vector<map_array> fast");
304     bench_my_outer_prod<ublas::sparse_matrix<T, ublas::row_major, ublas::map_array<std::size_t, T> >,
305                         ublas::sparse_vector<T, ublas::map_array<std::size_t, T> >, N> () (runs, fast_tag ());
306 #endif
307 
308 #ifdef USE_STD_MAP
309     header ("sparse_matrix<std::map>, sparse_vector<std::map> safe");
310     bench_my_outer_prod<ublas::sparse_matrix<T, ublas::row_major, std::map<std::size_t, T> >,
311                         ublas::sparse_vector<T, std::map<std::size_t, T> >, N> () (runs, safe_tag ());
312 
313     header ("sparse_matrix<std::map>, sparse_vector<std::map> fast");
314     bench_my_outer_prod<ublas::sparse_matrix<T, ublas::row_major, std::map<std::size_t, T> >,
315                         ublas::sparse_vector<T, std::map<std::size_t, T> >, N> () (runs, fast_tag ());
316 #endif
317 #endif
318 
319 #ifdef USE_COMPRESSED_MATRIX
320     header ("compressed_matrix, compressed_vector safe");
321     bench_my_outer_prod<ublas::compressed_matrix<T, ublas::row_major>,
322                         ublas::compressed_vector<T>, N> () (runs, safe_tag ());
323 
324     header ("compressed_matrix, compressed_vector fast");
325     bench_my_outer_prod<ublas::compressed_matrix<T, ublas::row_major>,
326                         ublas::compressed_vector<T>, N> () (runs, fast_tag ());
327 #endif
328 
329 #ifdef USE_COORDINATE_MATRIX
330     header ("coordinate_matrix, coordinate_vector safe");
331     bench_my_outer_prod<ublas::coordinate_matrix<T, ublas::row_major>,
332                         ublas::coordinate_vector<T>, N> () (runs, safe_tag ());
333 
334     header ("coordinate_matrix, coordinate_vector fast");
335     bench_my_outer_prod<ublas::coordinate_matrix<T, ublas::row_major>,
336                         ublas::coordinate_vector<T>, N> () (runs, fast_tag ());
337 #endif
338 
339 #ifdef USE_STD_VALARRAY
340     header ("std::valarray");
341     bench_cpp_outer_prod<std::valarray<T>, std::valarray<T>, N> () (runs);
342 #endif
343 
344     header ("prod (matrix, vector)");
345 
346     header ("C array");
347     bench_c_matrix_vector_prod<T, N> () (runs);
348 
349 #ifdef USE_SPARSE_MATRIX
350 #ifdef USE_MAP_ARRAY
351     header ("sparse_matrix<map_array>, sparse_vector<map_array> safe");
352     bench_my_matrix_vector_prod<ublas::sparse_matrix<T, ublas::row_major, ublas::map_array<std::size_t, T> >,
353                                 ublas::sparse_vector<T, ublas::map_array<std::size_t, T> >, N> () (runs, safe_tag ());
354 
355     header ("sparse_matrix<map_array>, sparse_vector<map_array> fast");
356     bench_my_matrix_vector_prod<ublas::sparse_matrix<T, ublas::row_major, ublas::map_array<std::size_t, T> >,
357                                 ublas::sparse_vector<T, ublas::map_array<std::size_t, T> >, N> () (runs, fast_tag ());
358 #endif
359 
360 #ifdef USE_STD_MAP
361     header ("sparse_matrix<std::map>, sparse_vector<std::map> safe");
362     bench_my_matrix_vector_prod<ublas::sparse_matrix<T, ublas::row_major, std::map<std::size_t, T> >,
363                                 ublas::sparse_vector<T, std::map<std::size_t, T> >, N> () (runs, safe_tag ());
364 
365     header ("sparse_matrix<std::map>, sparse_vector<std::map> fast");
366     bench_my_matrix_vector_prod<ublas::sparse_matrix<T, ublas::row_major, std::map<std::size_t, T> >,
367                                 ublas::sparse_vector<T, std::map<std::size_t, T> >, N> () (runs, fast_tag ());
368 #endif
369 #endif
370 
371 #ifdef USE_COMPRESSED_MATRIX
372     header ("compressed_matrix, compressed_vector safe");
373     bench_my_matrix_vector_prod<ublas::compressed_matrix<T, ublas::row_major>,
374                                 ublas::compressed_vector<T>, N> () (runs, safe_tag ());
375 
376     header ("compressed_matrix, compressed_vector fast");
377     bench_my_matrix_vector_prod<ublas::compressed_matrix<T, ublas::row_major>,
378                                 ublas::compressed_vector<T>, N> () (runs, fast_tag ());
379 #endif
380 
381 #ifdef USE_COORDINATE_MATRIX
382     header ("coordinate_matrix, coordinate_vector safe");
383     bench_my_matrix_vector_prod<ublas::coordinate_matrix<T, ublas::row_major>,
384                                 ublas::coordinate_vector<T>, N> () (runs, safe_tag ());
385 
386     header ("coordinate_matrix, coordinate_vector fast");
387     bench_my_matrix_vector_prod<ublas::coordinate_matrix<T, ublas::row_major>,
388                                 ublas::coordinate_vector<T>, N> () (runs, fast_tag ());
389 #endif
390 
391 #ifdef USE_STD_VALARRAY
392     header ("std::valarray");
393     bench_cpp_matrix_vector_prod<std::valarray<T>, std::valarray<T>, N> () (runs);
394 #endif
395 
396     header ("matrix + matrix");
397 
398     header ("C array");
399     bench_c_matrix_add<T, N> () (runs);
400 
401 #ifdef USE_SPARSE_MATRIX
402 #ifdef USE_MAP_ARRAY
403     header ("sparse_matrix<map_array> safe");
404     bench_my_matrix_add<ublas::sparse_matrix<T, ublas::row_major, ublas::map_array<std::size_t, T> >, N> () (runs, safe_tag ());
405 
406     header ("sparse_matrix<map_array> fast");
407     bench_my_matrix_add<ublas::sparse_matrix<T, ublas::row_major, ublas::map_array<std::size_t, T> >, N> () (runs, fast_tag ());
408 #endif
409 
410 #ifdef USE_STD_MAP
411     header ("sparse_matrix<std::map> safe");
412     bench_my_matrix_add<ublas::sparse_matrix<T, ublas::row_major, std::map<std::size_t, T> >, N> () (runs, safe_tag ());
413 
414     header ("sparse_matrix<std::map> fast");
415     bench_my_matrix_add<ublas::sparse_matrix<T, ublas::row_major, std::map<std::size_t, T> >, N> () (runs, fast_tag ());
416 #endif
417 #endif
418 
419 #ifdef USE_COMPRESSED_MATRIX
420     header ("compressed_matrix safe");
421     bench_my_matrix_add<ublas::compressed_matrix<T, ublas::row_major>, N> () (runs, safe_tag ());
422 
423     header ("compressed_matrix fast");
424     bench_my_matrix_add<ublas::compressed_matrix<T, ublas::row_major>, N> () (runs, fast_tag ());
425 #endif
426 
427 #ifdef USE_COORDINATE_MATRIX
428     header ("coordinate_matrix safe");
429     bench_my_matrix_add<ublas::coordinate_matrix<T, ublas::row_major>, N> () (runs, safe_tag ());
430 
431     header ("coordinate_matrix fast");
432     bench_my_matrix_add<ublas::coordinate_matrix<T, ublas::row_major>, N> () (runs, fast_tag ());
433 #endif
434 
435 #ifdef USE_STD_VALARRAY
436     header ("std::valarray");
437     bench_cpp_matrix_add<std::valarray<T>, N> () (runs);
438 #endif
439 }
440 
441 #ifdef USE_FLOAT
442 template struct bench_2<float, 3>;
443 template struct bench_2<float, 10>;
444 template struct bench_2<float, 30>;
445 template struct bench_2<float, 100>;
446 #endif
447 
448 #ifdef USE_DOUBLE
449 template struct bench_2<double, 3>;
450 template struct bench_2<double, 10>;
451 template struct bench_2<double, 30>;
452 template struct bench_2<double, 100>;
453 #endif
454 
455 #ifdef USE_STD_COMPLEX
456 #ifdef USE_FLOAT
457 template struct bench_2<std::complex<float>, 3>;
458 template struct bench_2<std::complex<float>, 10>;
459 template struct bench_2<std::complex<float>, 30>;
460 template struct bench_2<std::complex<float>, 100>;
461 #endif
462 
463 #ifdef USE_DOUBLE
464 template struct bench_2<std::complex<double>, 3>;
465 template struct bench_2<std::complex<double>, 10>;
466 template struct bench_2<std::complex<double>, 30>;
467 template struct bench_2<std::complex<double>, 100>;
468 #endif
469 #endif
470