# Ublas :

From: Karl Meerbergen (Karl.Meerbergen_at_[hidden])
Date: 2007-03-02 09:12:05

Sourabh wrote:

>Hi,
>Can anyone please suggest me some optimizations in the loop to improve
>the overall speed of the following program ?
>
>The steps I am doing are:
>1. Solve for Y the equation: (I-X) Y = B, where X is a lower triangular
>sparse matrix.
>
>

It should be possible to do this with the triangular solves from ublas
but I am unsure whether they efficiently handle sparse matrices. Gunter,

I suggest you store I-X in a sparse matrix and use a triangular solve
routine.

>2. Multiply Y with lower triangular sparse matrix Z, A = Y*Z.
>
>

I do not know which is the best routine for that: not prod() I guess.

>3. Add A to another vector.
>
>

use other_vector.plus_assign(a)

>The sizes in the above steps are of order of 100.
>
>The above mentioned sequence is repeated by 2e9 times. So it is taking too
>much time.
>I am sure if I optimize the inner-most loop, it can be made fast.
>What optimizations can be made in the loop ?
>
>

Are the triangular sparse matrices very sparse, or almost dense? If they
are almost dense, use the BLAS bindings for dense matrices.

Karl

>
>#include <boost/numeric/ublas/matrix.hpp>
>#include <boost/numeric/ublas/triangular.hpp>
>#include <boost/numeric/ublas/io.hpp>
>#include <boost/numeric/ublas/matrix_sparse.hpp>
>#include <boost/numeric/ublas/matrix_proxy.hpp>
>using namespace boost::numeric::ublas;
>
>int main (int argc, char* argv[])
>{
> unsigned int samples = atoi (argv[1]);
> unsigned int numpaths = atoi (argv[2]);
> unsigned short gatesperpath = atoi (argv[3]);
>
> for (unsigned int i = 0; i < samples; ++i) {
> for (unsigned int j = 0; j < numpaths ; ++ j) {
> mapped_matrix<double> lembda (gatesperpath, gatesperpath);
> identity_matrix<double> ident (gatesperpath);
> slopes = solve (lembda, omegas, lower_tag ());
> vector <double> delays (gatesperpath);
> delays = prod (lembda, slopes);
> delays += omegas;
> }
> }
>
>}
>
>