Boost logo

Ublas :

From: David M Garza (David.M.Garza_at_[hidden])
Date: 2006-02-15 12:22:03

I'm working on a code where I'm building up the matrix from many small
dense blocks. It's not a finite element code, but I would think it's
very similar in its pattern of building up a large sparse matrix from
small dense ones. I was thinking it would be nicer to write something like

   subrange(A,i0,i0 + m,j0,j0 + n) += b

rather than

   for(int i = 0;i < m;i++) {
     for(int j = 0;j < n;j++) {
       A.append_element(i0 + i,j0 + j,b(i,j))

Since the coordinate matrix supports summation of repeated elements I
thought the += and -= operators might be written to take adavntage of
this fact and just use append_element().

I got as far as the matrix_assign<scalar_plus_assign> in the source and
didn't follow it further to figure out if that line would perhaps lead
to an optimized routine. Would replacing the matrix_assign line with
the nested loop above be the right place to introduce the optimization?

Once I have the sparse matrix built I hand it over to an optimizer which
expects the matrix to be in compressed sparse column form. My question
above came from the larger issue of which is faster--assembling a
coordinate_matrix and converting to a compressed_matrix at the end, or
assembling using the compressed_matrix.

I've seen your pages on sparse matrix comparisons, and filling a
compressed matrix with a preexisting structure is my current plan since
I have to give the structure to the optimizer when I first call it. I'm
dealing with smaller problems than you tested with--at the large end
I'll have a matrix which is something like 30000x40000 with around
300000 nonzero entries--and until I build and test the code I won't have
a good feel for where the speed issues are. I could easily be worrying
about something which isn't a problem.


Gunter Winkler wrote:
> On Tuesday 14 February 2006 17:30, David M Garza wrote:
>>Are the overloaded += and -= operators for coordinate_matrix set up to
>>simply append the elements to the end of the matrix arrays?
> Looking at matrix_sparse.hpp around line 4232 shows that the += (plus_assign),
> -= (minus_assign) are not optimized in any way. They use the
> matrix_assign<scalar_plus_assign> (*this, ae);
> algorithm.
> Personally, I never used these operators. What kind of application do you have
> in mind?
> mfg
> Gunter
> _______________________________________________
> ublas mailing list
> ublas_at_[hidden]