<em>guwi17_at_[hidden] wrote:<br><br></em>>BTW. I would not expect a big difference between a banded matrix with <br>>compile time known bandwidth compared to the current implementation. I <br>>still wait for someone who proves the opposite ;-)<br>><br>>Do you already have a benchmark to measure the speed of banded_matrix? <br><br>I am doing a very simple test with the (undocumented?) diagonal_matrix in banded.hpp<br><br>A matrix product of two diagonal_matrix should return another diagonal_matrix and should have linear complexity O(size). The same goes for an element-wise vector product. Storage should not be an issue, as the data layout of the diagonal_matrix should match that of the vector. <br><br>The test performs one of both, and time measurements for a size of 1E7 double elements shows that the vector product is 10 times faster. The test also uses a matrix_vector_slice to get a view of the diagonal_matrix as a vector and does an element-wise product. This is much faster that the matrix product, but vector product is still twice as fast. <br><br>IMO compile time knowledge of the bandwith of a banded_matrix *IS* important, specially for large sizes and small bandwidth (being a large diagonal_matrix the extreme case).<br><br>I suppose that implementing a banded_matrix_c with compile time bandwith is not trivial, at least not if optimal traversal of elements is required... <br><br><br><br>test code:<br><br><br>#include <boost/numeric/ublas/matrix.hpp><br>#include <boost/numeric/ublas/matrix_proxy.hpp><br>#include <boost/numeric/ublas/banded.hpp><br><br>#ifdef NDEBUG<br> const int SIZE = 10000000; <br>#else<br> const int SIZE = 10;<br>#endif<br><br>void diagonal()<br>{<br> using namespace boost::numeric::ublas;<br> diagonal_matrix<double> d1 ( SIZE );<br> diagonal_matrix<double> d2 ( SIZE ); <br> <br> // how to do this with diagonal_adaptor<scalar_matrix> ?<br> std::fill(d1.data().begin(), d1.data().end(), 2.0);<br> std::fill(d2.data().begin(), d2.data().end(), 10.0);<br> <br> diagonal_matrix<double> d3 = prod(d1,d2); <br>}<br><br>void diagonal_as_vector()<br>{<br> using namespace boost::numeric::ublas;<br> diagonal_matrix<double> d1 ( SIZE );<br> diagonal_matrix<double> d2 ( SIZE );<br> <br> // how to do this with diagonal_adaptor<scalar_matrix> ? <br> std::fill(d1.data().begin(), d1.data().end(), 2.0);<br> std::fill(d2.data().begin(), d2.data().end(), 10.0);<br><br> matrix_vector_slice< diagonal_matrix<double> > dx1(d1, slice(0,1,SIZE), slice(0,1,SIZE)); <br> matrix_vector_slice< diagonal_matrix<double> > dx2(d2, slice(0,1,SIZE), slice(0,1,SIZE));<br> <br> // any vector_expression to diagonal_matrix ?<br> vector<double> d3 = element_prod(dx1, dx2); <br>}<br><br>void vectors()<br>{<br> using namespace boost::numeric::ublas;<br> vector<double> v1(scalar_vector<double>(SIZE, 2.0));<br> vector<double> v2(scalar_vector<double>(SIZE, 10.0)); <br> vector<double> v3 = element_prod(v1, v2);<br>}<br><br>int main(int argc, char* argv[])<br>{<br> if (argc != 2) return 1;<br> if (*argv[1] == 'd') diagonal();<br> else if (*argv[1] == 'v') vectors(); <br> else if (*argv[1] == 'x') diagonal_as_vector();<br> return 0;<br>}<br><br>