 # Ublas :

From: Regis Behmo (regisb_at_[hidden])
Date: 2007-09-03 12:10:08

On 9/3/07, Paul C. Leopardi <paul.leopardi_at_[hidden]> wrote:
> On Tue, 4 Sep 2007, Regis Behmo wrote:
> > Hello everyone again,
> > How would you iterate over the i-th row or the j-th column of a
> > mapped_matrix or a compressed_matrix? I got the feeling it's not possible
> > to do it in a proper (i.e: optimal) way.
> >
> > In the case of a simple matrix, I would do something like this:
> >
> > iterator1 it1 = this->begin1();
> > for( int x = 0; x < i; x++ )
> > it1++;
> > iterator it2 = it1.begin();
> >
> >
> > Régis
>
> Hi Regis,
> What is i?

i is the index of the row I want to iterate over.

> Why does your loop not check it1 != this->end1()?

Because I suppose i is inferior than the number of rows of my matrix

> What will you be using x for?

The only thing x is used for is to make sure I reach the i-th row in
the for loop.

I will give another (clearer) example:

int main()
{
ublas::matrix M(10,10); // I create a matrix with 10 rows and 10 columns
int i = 4; // I want to print the elements from the 4th row

iterator1 it1 = M.begin1();
for( int x = 0; x < i; x++ )
it1++;
iterator it2 = it1.begin();
while( it2 != it1.end() )
{
cout << *it2 << endl;
it2++;
}

}

This example prints the elements of the i-th row. Now if I want to do
this for a sparse matrix, things get more complex, don't they? For
instance, if some rows of my matrix M are empty, the first for loop is
going to reach a row that is not the i-th one. I would have to check
the index of it1 in order to do that. I wonder if there is any more
appropriate way do this.

> Following is an example from GluCat which apparently works for both dense and
> compressed matrices.
> Best, Paul
>
> /// Kronecker tensor product of matrices - as per Matlab kron
> template< typename LHS_T, typename RHS_T >
> const
> RHS_T
> kron(const LHS_T& lhs, const RHS_T& rhs)
> {
> typedef RHS_T matrix_t;
> typedef typename matrix_t::size_type matrix_index_t;
> const matrix_index_t rhs_s1 = rhs.size1();
> const matrix_index_t rhs_s2 = rhs.size2();
> matrix_t result(lhs.size1()*rhs_s1, lhs.size2()*rhs_s2);
> result.clear();
> typedef typename matrix_t::value_type Scalar_T;
> typedef typename LHS_T::const_iterator1 lhs_const_iterator1;
> typedef typename LHS_T::const_iterator2 lhs_const_iterator2;
> typedef typename RHS_T::const_iterator1 rhs_const_iterator1;
> typedef typename RHS_T::const_iterator2 rhs_const_iterator2;
> for (lhs_const_iterator1
> lhs_it1 = lhs.begin1();
> lhs_it1 != lhs.end1();
> ++lhs_it1)
> for (lhs_const_iterator2
> lhs_it2 = lhs_it1.begin();
> lhs_it2 != lhs_it1.end();
> ++lhs_it2)
> {
> const matrix_index_t start1 = rhs_s1 * lhs_it2.index1();
> const matrix_index_t start2 = rhs_s2 * lhs_it2.index2();
> const Scalar_T& lhs_val = (*lhs_it2);
> for (rhs_const_iterator1
> rhs_it1 = rhs.begin1();
> rhs_it1 != rhs.end1();
> ++rhs_it1)
> for (rhs_const_iterator2
> rhs_it2 = rhs_it1.begin();
> rhs_it2 != rhs_it1.end();
> ++rhs_it2)
> result(start1 + rhs_it2.index1(), start2 + rhs_it2.index2()) =
> lhs_val * *rhs_it2;
> }
> return result;
> }
>

In this example there is no iterator on result, there are only
iterators on all lines and columns of lhs and rhs.

Régis