boost compressed matrix basics

Question

I am confused on how the boost::compressed_matrix works. Suppose I declare the compressed_matrix like this:

boost::numeric::ublas::compressed_matrix<double> T(1000, 1000, 3*1000);

This allocates space for 3*1000 elements in a 1000x1000 matrix. Now how do I give it the locations which are the non-zero elements? When and how are the non-zero elements set? Is it each time I assign an element in the matrix, e.g. B(4,4)=4, it would mark that element as non-zero?

I would really appreciate if you could help me learn this using an example if possible. Some insight into the internal implementation would be great. I want to make sure I don't write programs that are sub-optimal by guess work.

thank you!

Answer 1

you can either use (i,j) operator to create nonzero element implicitly or use insert_element function to insert element explicitly.

Best place is actually look inside implementation:

http://www.tena-sda.org/doc/5.2.2/boost/d2/db7/matrix__sparse_8hpp-source.html#l02761

true_reference insert_element (size_type i, size_type j, const_reference t)

Inserts the value t at the j-th element of the i-th row. Duplicates elements are not allowed.

---

void erase_element (size_type i, size_type j)

Erases the value at the j-th element of the i-th row.

Answer 2

compressed matrix has an underlying linear container (unbounded_array by default, but you can make it bounded_array or std::vector if you want), which contains all non-zero elements of the matrix, in row-major (by default) order. That means that whenever you write a new non-zero element to compressed matrix, it is inserted into that underlying array. If you're not filling the matrix in (row-major) order, every insert will be O(n). When you're changing an existing non-zero element, it is simply changed in the underlying array.

Here's a simple test to see what the underlying structure looks like:

#include <boost/numeric/ublas/matrix_sparse.hpp>
#include <boost/numeric/ublas/storage.hpp>
namespace ublas = boost::numeric::ublas;
void show_array(const ublas::unbounded_array<double>& a)
{
    for(size_t i=0; i<a.size(); ++i)
            std::cout << a[i] << ' ';
    std::cout << '\n';
}
int main()
{
    ublas::compressed_matrix<double> m (10, 10, 3 * 10);
    m(0, 5) = 1; // underlying array is {1, 0, 0, 0, ...}
    show_array(m.value_data());
    m(0, 6) = 2; // underlying array is {1, 2, 0, 0, ...}
    show_array(m.value_data());
    m(0, 4) = 3;  // underlying array is {3, 1, 2, 0, ...}
    show_array(m.value_data());
    m(0, 4) = 7;  // underlying array is {7, 1, 2, 0, ...}
    show_array(m.value_data());
}