|
|
Ublas : |
From: Gunter Winkler (guwi17_at_[hidden])
Date: 2006-05-24 03:08:51
On Tuesday 23 May 2006 20:39, eeitan_at_[hidden] wrote:
> thanks for your fast answear.
> i am tring to do "cartesian product" of matrices:
> matA =[2,2]((a11,a12),
> (a21,a22))
> matB =[2,2]((b11,b12),
> (b21,b22))
> so the cartesian product of them is
> matC = [4,4]((a11*matB,a12*matB),
> (a21*matB,a22*matB))
> so i made this function:
>
> // making a cartesian spaces product between matrices
> template< class M1,class M2,class M3 >
> void cartesianSpacesProduct( M1 & target ,const M2 & left ,const M3 &
> right) {// not deleting the target at this function !!!!!!!!
> assert(target.size1()==left.size1()*right.size1() &&
> target.size2()==left.size2()*right.size2() &&
> target.size1()==target.size2() &&
> left.size1()==left.size2());
> for(typename M2::const_iterator1 it1(left.begin1());it1 != left.end1()
> ; ++it1){
> for (typename M2::const_iterator2 it2(it1.begin());it2 !=
> it1.end(); ++it2){
> ublas::matrix_range<M1 > mr (target,
> ublas::range
> (it2.index1()*right.size1()
>
> ,(it2.index1()+1)*right.size1()),
> ublas::range
> (it2.index2()*right.size2()
>
> ,(it2.index2()+1)*right.size2()));
> mr+=(*it2)*right;
> }
> }
> };
> this function work very slow for me what you recommend me to do to boost it
> up ?
You can possibly speed up the assignment if you use
noalias( mr ) += (*it2) * right;
For such expressions usually a temporary matrix is computed, because there is
no easy way to decide whether there is aliasing or not.
For dense matrices this should be quite fast, so I suspect that you multiply
sparse matrices. If the target is sparse, then selecting a submatrix is very
slow, because the non-const operator(i,j) or the find() method is called many
times. IMHO it is hard to tell where the most time is lost without running a
profiler. Maybe you have a look at the different axpy_prods in operation.hpp
and write a specialized Kronecker_prod.
For optimal performance you have to combine two things:
* Fill the lhs matrix in order (e.g. with push_back)
* Use iterators at rhs
Did you try to compute the lhs column by column or row by row?
mfg
Gunter