$include_dir="/home/hyper-archives/ublas/include"; include("$include_dir/msg-header.inc") ?>
Subject: Re: [ublas] matrix_range for matrix
From: Nasos Iliopoulos (nasos_i_at_[hidden])
Date: 2013-07-05 16:44:38
Petro,
here is a compilable and working version. GCC 4.7. had many issues with 
the previous code. What compiler are you using? You can probably carve 
this to fit your needs but I would suggest to rethink your approach.
Let me know if it works.
Best,
-Nasos
#include <boost/numeric/ublas/storage.hpp>
#include <boost/numeric/ublas/matrix.hpp>
#include <boost/numeric/ublas/matrix_proxy.hpp>
#include <boost/numeric/ublas/triangular.hpp>
using namespace boost::numeric;
/*    Really a+ boost general matrix_range to
    a properly-lda-sized boost matrix it
    holds internally.
*/
template < typename _T, typename _L = ublas::column_major, typename _A = 
ublas::unbounded_array<_T> >
class MatrixGeneralT
    : public ublas::matrix_expression< MatrixGeneralT< _T, _L, _A>  >
{
    size_t ld_ ;
    size_t
        determineLDA (
            const size_t nRows
            , const bool mklLeadingDimension = true
    ) {
        static const size_t MKL_LD_DIVSOR_IN_BYTES = 32 ;
        static const size_t MKL_LD_DIVISOR_DISALLOWED = 2028 ;
        if ( ! mklLeadingDimension )
            return ( ld_ = nRows ) ;
        ld_ = nRows ;
    // the nbr of bytes per lda should be divisible by 
MKL::LDA_DIVISOR_IN_BYTES.
        size_t bytesNbrRows =
            nRows * sizeof( _T ) ;
        const ldiv_t qr =
            ldiv( long(bytesNbrRows), long( MKL_LD_DIVSOR_IN_BYTES ) ) ;
    //if not, resize lda so that this is true.
        if ( qr.rem != 0 )
            bytesNbrRows =
                ( qr.quot + 1 ) * MKL_LD_DIVSOR_IN_BYTES ;
    // but avoid the bad dimension ( e.g. 2048)
        while ( bytesNbrRows % MKL_LD_DIVISOR_DISALLOWED == 0 )
            bytesNbrRows += MKL_LD_DIVSOR_IN_BYTES ;
        return
            ( ld_ = bytesNbrRows / sizeof( _T ) ) ;
    }
public:
    typedef
        MatrixGeneralT                self_type ;
    typedef typename
        ublas::matrix< _T, _L, _A >                M ;
    typedef
        M internal_matrix_type ;
    typedef
        const M const_internal_matrix_type ;
    typedef typename
        ublas::matrix_range<
            internal_matrix_type
        > matrix_range_type ;
    typedef typename
        ublas::matrix_range<
            const_internal_matrix_type
        > const_matrix_range_type ;
    typedef     M                        matrix_type ;
    typedef typename    M::size_type            size_type ;
    typedef typename    M::difference_type        difference_type ;
    typedef typename    M::value_type            value_type ;
    typedef typename    M::const_reference        const_reference ;
    typedef typename std::conditional<
                        std::is_const<M>::value,
                        typename M::const_reference,
                        typename M::reference
    >::type reference ;
    typedef typename std::conditional<
                        std::is_const<M>::value ,
                        typename M::const_closure_type,
                        typename M::closure_type
                    >::type matrix_closure_type ;
    typedef ublas::basic_range<
                        size_type
                        , difference_type
            > range_type ;
    typedef
        const_matrix_range_type const_closure_type;
    typedef
        matrix_range_type closure_type;
    typedef typename
        ublas::storage_restrict_traits<
                typename M::storage_category,
                ublas::dense_proxy_tag
        >::storage_category storage_category ;
    typedef typename
        M::orientation_category orientation_category;
    typedef typename
        matrix_range_type::iterator1 iterator1 ;
    typedef typename
        matrix_range_type::iterator2 iterator2 ;
    typedef typename
        const_matrix_range_type::iterator1 const_iterator1 ;
    typedef typename
        const_matrix_range_type::iterator2 const_iterator2 ;
private :
    size_type    size1_ ;
    size_type    size2_ ;
    internal_matrix_type    matrix_ ;
public :
    MatrixGeneralT( bool optimalLDA = true )
        :size1_(0)
        ,size2_(0)
        , matrix_(0,0)
    {}
    MatrixGeneralT(
        const size_t nbrRows,
        const size_t nbrCols = 0
    )
        :size1_( nbrRows )
        ,size2_( nbrCols )
        , matrix_( determineLDA( nbrRows ), nbrCols )
    {}
       MatrixGeneralT(
        const size_t nbrRows
    )
        :size1_( nbrRows )
        ,size2_( nbrRows )
        , matrix_( determineLDA( nbrRows ), nbrRows )
    {}
    virtual ~MatrixGeneralT()
    {}
// auto casting
    operator matrix_range_type () {
        return  matrix_range_type( matrix_
                , ublas::range( 0, size1_ )
                , ublas::range( 0, size2_ )
            ) ;
    }
    operator const_matrix_range_type () const {
        return const_matrix_range_type ( matrix_, ublas::range( 0, 
size1_ ), ublas::range( 0, size2_ ) ) ;
    }
// matrix_expression assignment :
    template <typename _E >
    self_type &
    operator = ( typename ublas::matrix_expression<_E> const & me ) {
        _E const & e = me() ;
        matrix_.resize( determineLDA( e.size1(), true ), e.size2() ) ;
        size1_ = e.size1() ;
        size2_ = e.size2() ;
        matrix_range_type
        (
            matrix_
            , ublas::range( 0, size1_ )
            , ublas::range( 0, size2_ )
        )
        =
        me() ;
        return
            *this ;
    }
    value_type & operator() ( const size_t i, const size_t j ) {
        return
            matrix_.data()[ i + ld_ * j ] ;
    }
    value_type const & operator() ( const size_t i, const size_t j ) const {
        return
            matrix_.data()[ i + ld_ * j ] ;
    }
} ;
int main() {
    typedef MatrixGeneralT< double  >    matrix ;
    matrix m(3,3) ;
    ublas::vector< double > x( 3 ) ;
    x[0] = 1; x[1] = 2 ; x[2] = 3;
    m =  ublas::outer_prod( x,x ) ;
    matrix m1( 3, 4 ) ;
    matrix m2( 4, 5 ) ;
    matrix m3( 3, 5 ) ;
    for ( size_t i = 0; i != 3 ; ++i )
        for ( size_t j = 0; j != 4 ; ++j )
            m1(i,j) = (i+1)  + (j+1)*10 ;
    for ( size_t i = 0; i != 4 ; ++i )
        for ( size_t j = 0; j != 5 ; ++j )
            m2(i,j) = (i+1)  + (j+1)*10 ;
    m3 = ublas::prod( m1, m2 ) ;
    return
        1;
}
On 07/05/2013 03:40 PM, Petros wrote:
> Naso,
>
> This is the skeleton of the code, in case someone needs a starting point.
> Please, rest assured that there are shortcomings and pitfalls but vary 
> likely
> a better starting point than nothing:
>
> //boost
> #include <boost/numeric/ublas/storage.hpp>
> #include <boost/numeric/ublas/matrix.hpp>
> #include <boost/numeric/ublas/matrix_proxy.hpp>
> //#include <boost/numeric/ublas/triangular.hpp>
> namespace ublas =  boost::numeric::ublas;
>
> /*    Really a+ boost general matrix_range to
>    a properly-lda-sized boost matrix it
>    holds internally.
> */
> template < typename _T, typename _L = ublas::column_major, typename _A 
> = ublas::unbounded_array<_T> >
> class MatrixGeneralT
>    : public ublas::matrix_expression< MatrixGeneralT< _T, _L, _A>  >
> {
>    size_t ld_ ;
>    size_t
>        determineLDA (
>            const size_t nRows
>            , const bool mklLeadingDimension = true
>    ) {
>        static const size_t MKL_LD_DIVSOR_IN_BYTES = 32 ;
>        static const size_t MKL_LD_DIVISOR_DISALLOWED = 2028 ;
>        if ( ! mklLeadingDimension )
>            return ( ld_ = nRows ) ;
>        ld_ = nRows ;
>    // the nbr of bytes per lda should be divisible by 
> MKL::LDA_DIVISOR_IN_BYTES.
>        size_t bytesNbrRows =
>            nRows * sizeof( _T ) ;
>        const ldiv_t qr =
>            div( long(bytesNbrRows), long( MKL_LD_DIVSOR_IN_BYTES ) ) ;
>    //if not, resize lda so that this is true.
>        if ( qr.rem != 0 )
>            bytesNbrRows =
>                ( qr.quot + 1 ) * MKL_LD_DIVSOR_IN_BYTES ;
>    // but avoid the bad dimension ( e.g. 2048)
>        while ( bytesNbrRows % MKL_LD_DIVISOR_DISALLOWED == 0 )
>            bytesNbrRows += MKL_LD_DIVSOR_IN_BYTES ;
>        return
>            ( ld_ = bytesNbrRows / sizeof( _T ) ) ;
>    }
>
> public:
>    typedef typename
>        MatrixGeneralT< _T, _L, _A>                self_type ;
>    typedef typename
>        ublas::matrix< _T, _L, _A >                M ;
>    typedef typename
>        M internal_matrix_type ;
>    typedef typename
>        const M const_internal_matrix_type ;
>    typedef typename
>        ublas::matrix_range<
>            internal_matrix_type
>        > matrix_range_type ;
>    typedef typename
>        ublas::matrix_range<
>            const_internal_matrix_type
>        > const_matrix_range_type ;
>
>    typedef typename    M                        matrix_type ;
>    typedef typename    M::size_type            size_type ;
>    typedef typename    M::difference_type        difference_type ;
>    typedef typename    M::value_type            value_type ;
>    typedef typename    M::const_reference        const_reference ;
>
>    typedef typename std::conditional<
>                        std::is_const<M>::value,
>                        typename M::const_reference,
>                        typename M::reference
>    >::type reference ;
>    typedef typename std::conditional<
>                        std::is_const<M>::value ,
>                        typename M::const_closure_type,
>                        typename M::closure_type
>                    >::type matrix_closure_type ;
>    typedef ublas::basic_range<
>                        size_type
>                        , difference_type
>            > range_type ;
>    typedef typename
>        const_matrix_range_type const_closure_type;
>    typedef typename
>        matrix_range_type closure_type;
>    typedef typename
>        ublas::storage_restrict_traits<
>                typename M::storage_category,
>                ublas::dense_proxy_tag
>        >::storage_category storage_category ;
>    typedef typename
>        M::orientation_category orientation_category;
>
>    typedef typename
>        matrix_range_type::iterator1 iterator1 ;
>    typedef typename
>        matrix_range_type::iterator2 iterator2 ;
>    typedef typename
>        const_matrix_range_type::iterator1 const_iterator1 ;
>    typedef typename
>        const_matrix_range_type::iterator2 const_iterator2 ;
>
> private :
>    typename internal_matrix_type    matrix_ ;
>    typename size_type    size1_ ;
>    typename size_type    size2_ ;
> public :
>
>    MatrixGeneralT( bool optimalLDA = true )
>        :size1_(0)
>        ,size2_(0)
>        , matrix_(0,0)
>    {}
>    MatrixGeneralT(
>        const size_t nbrRows,
>        const size_t nbrCols = 0
>    )
>        :size1_( nbrRows )
>        ,size2_( nbrCols == 0 ? size1_ : nbrCols )
>        , matrix_( determineLDA( nbrRows ), nbrCols )
>    {}
>    virtual ~MatrixGeneralT()
>    {}
> // auto casting
>    operator matrix_range_type () {
>        return
>            matrix_range_type
>            (
>                matrix_
>                , ublas::range( 0, size1_ )
>                , ublas::range( 0, size2_ )
>            ) ;
>    }
>    operator const_matrix_range_type () const {
>        return
>            const_matrix_range_type
>            (
>                matrix_
>                , ublas::range( 0, size1_ )
>                , ublas::range( 0, size2_ )
>            ) ;
>    }
>
> // matrix_expression assignment :
>    template <typename _E >
>    self_type &
>    operator = ( typename ublas::matrix_expression<_E> const & me ) {
>        _E const & e = me() ;
>        matrix_.resize( determineLDA( e.size1(), true ), e.size2() ) ;
>        size1_ = e.size1() ;
>        size2_ = e.size2() ;
>        matrix_range_type
>        (
>            matrix_
>            , ublas::range( 0, size1_ )
>            , ublas::range( 0, size2_ )
>        )
>        =
>        me() ;
>        return
>            *this ;
>    }
>
>    value_type & operator() ( const size_t i, const size_t j ) {
>        return
>            matrix_.data()[ i + ld_ * j ] ;
>    }
>    value_type const & operator() ( const size_t i, const size_t j ) 
> const {
>        return
>            matrix_.data()[ i + ld_ * j ] ;
>    }
>
> } ;
>
> int main() {
>
>    typedef MatrixGeneralT< double  >    matrix ;
>
>    matrix m(3,3) ;
>    ublas::vector< double > x( 3 ) ;
>    x[0] = 1; x[1] = 2 ; x[2] = 3;
>    m =  ublas::outer_prod( x,x ) ;
>
>    matrix m1( 3, 4 ) ;
>    matrix m2( 4, 5 ) ;
>    matrix m3( 3, 5 ) ;
>
>    for ( size_t i = 0; i != 3 ; ++i )
>        for ( size_t j = 0; j != 4 ; ++j )
>            m1(i,j) = (i+1)  + (j+1)*10 ;
>    for ( size_t i = 0; i != 4 ; ++i )
>        for ( size_t j = 0; j != 5 ; ++j )
>            m2(i,j) = (i+1)  + (j+1)*10 ;
>    m3 = ublas::prod( m1, m2 ) ;
>
>    return
>        1;
> }
>
> -----Original Message----- From: Nasos Iliopoulos
> Sent: Friday, July 05, 2013 1:15 PM
> To: ublas_at_[hidden]
> Subject: Re: [ublas] matrix_range for matrix
>
> Two things:
> As part of GSoC Sathyam is building vector/Matrix view classes, along
> with aligned allocators, that might help you get a much more clean
> solution when interfacing with external code. His vector code and
> aligned allocators are  already in good shape and I expect we will get
> them into uBlas soon (as in 3 months).
>
> Additionally if you could provide some compilable code it would help
> debug the problem.
>
> Best,
> -Nasos
>
> On 07/05/2013 01:04 PM, pmamales_at_[hidden] wrote:
>> Naso,
>> Thx for getting back to me.
>>
>> A clarification : the matrix_range'ness of my matrix is coming into 
>> play only through the auto cast operator!
>> It is really smoke and mirrors. Tried to "pipe" my class into ublas 
>> by the typedefs and by deriving from
>> matrix_expression.
>>
>> The auto cast operation was also failing me for the assignment of a 
>> ublas::matrix_expression.
>> So, I said, fine, I will "view" the general matrix (matrix_) as a 
>> range and perform the assignment
>> there. What is returned is the object itself, not the "view".
>> Now why the assignment was failing me, I do not know either ;-)). But 
>> assignments are more tricky
>> and I thought of the ublas::prod test.
>>
>> All this is a hackery, unfortunately. A clean solution, i.e. one that 
>> would provide with a
>> leading dimension, decoupling the (fortran-layout matrices) leading 
>> dimension from size1() is not present.
>> Very likely, ublas designers wanted this behavior to be taken care by 
>> matrix_range, which is the more
>> appropriate c++ way.
>>
>> However, LAPACK has these "oddities" to account for the FORTRAN 
>> limitations ( at least back then ).
>> You see, I use ublas structures as containers, for lapack operations 
>> - speed is essential in what I do.
>> Occasionally, and for all kinds of test calculations ( or even in the 
>> private area of my clases,
>> for at least temporary purposes) I would like my struct's to allow 
>> for the ublas operations.
>>
>> I am stuck with this for a couple of days now, and try to find a 
>> clean way out.
>> The auto cast operation seemed like a good candidate. However, I get 
>> problems that I did not anticipate
>> and do not understand..
>> If you or anyone else could point to a better "take" I "am all ears" ;-)
>>
>> Thank you for your help
>> Petros
>>
>> ps: Very likely, everybody understands that, I forgot to include the 
>> ublas headers and the namespace alias,
>> in the code snippet I sent. Apologies..
>>
>> ---- Nasos Iliopoulos <nasos_i_at_[hidden]> wrote:
>>> Petro, Excuse me, but I hit send before finishing the e-mail:
>>>
>>> What are you trying to achieve here?:
>>>
>>>       template <typename _me>
>>>       matrix_ge_self_type  &
>>>       operator  = (_me const & me )
>>>       {
>>>               matrix_range_type
>>>               (
>>>                   matrix_
>>>                   , ublas::range( 0, size1_ )
>>>                   , ublas::range( 0, size2_ )
>>>               )  = me ;
>>>           return
>>>               *this ;
>>>       }
>>>
>>> it looks that you are creating a temporary that dies before ever being
>>> assigned to the encapsulated matrix_renge
>>>
>>> If this is the case than the matrix_range will be null, and hence 
>>> trying
>>> to write to it will probably crash everything.
>>>
>>> -Nasos
>>>
>>>
>>> On 07/05/2013 10:51 AM, pmamales_at_[hidden] wrote:
>>>> In order to accommodate for special memory layout for FORTRAN-type 
>>>> matrices I need to use a "matrix" which is really a 
>>>> ublas::matrix_range on a
>>>> ublas::matrix ( memory has to be aligned on 32-byte memory and 
>>>> actually every column of the matrix has to use this kind of 
>>>> boundary ).
>>>> In order to do this, I am defining my own ublas-conformant type 
>>>> (called MatrixGeneralT and declared as a ublas::matrix_expression) 
>>>> which has the same template arguments as the ublas::matrix and defines
>>>> the type-"interface" of a matrix range. To make things "happen" I 
>>>> create an automatic conversion from my class to the 
>>>> ublas::matrix_range.
>>>> However, when I try, as a first test, to calculate the product of 2 
>>>> matrices, the matrix_matrix_binary constructor fails to convert my 
>>>> class to matrix_range !
>>>> Unless, I have done something that c++ forbids ( which I don't 
>>>> think, since tried the "trick" with toy classes ) there must be 
>>>> something in ublas preventing this from
>>>> working.
>>>> Any ideas/input will be greatly appreciated.
>>>> TIA for your help,
>>>> Petros
>>>> ps: I use intel compiler 12.1 on windows 7 64bit plugged into 
>>>> msvc2010 ( hence the stl "flavor" ) and C++0x.
>>>>
>>>> I supply here a pseudo code for anyone interested in it :
>>>> #include <boost/numeric/ublas/matrix_expression.hpp>
>>>> template < typename _T, typename _L = ublas::column_major, typename 
>>>> _A = ublas::unbounded_array<_T> >
>>>> class MatrixGeneralT
>>>> :public ublas::matrix_expression< MatrixGeneralT<_T, _L, _A> >
>>>> {
>>>> size_t ld_ ;
>>>>
>>>> size_t
>>>> determineLDA (
>>>> const size_t nRows
>>>> , const bool mklLeadingDimension
>>>> ) {
>>>> if ( ! mklLeadingDimension )
>>>> return ( ld_ = nRows ) ;
>>>> ld_ = nRows ;
>>>> // the nbr of bytes per lda should be divisible by 
>>>> MKL::LDA_DIVISOR_IN_BYTES.
>>>> size_t bytesNbrRows =
>>>> nRows * sizeof( _T ) ;
>>>> const ldiv_t qr =
>>>> div( long(bytesNbrRows), long( 32/*MKL::LDA_DIVISOR_IN_BYTES*/ ) ) ;
>>>> //if not, resize lda so that this is true.
>>>> if ( qr.rem != 0 )
>>>> bytesNbrRows =
>>>> ( qr.quot + 1 ) * 32 /*MKL::LDA_DIVISOR_IN_BYTES */;
>>>> // but avoid the bad dimension ( e.g. 2048)
>>>> while ( bytesNbrRows % MKL::LDA_DIVISOR_IN_BYTES_FORBITTEN == 0 )
>>>> bytesNbrRows += MKL::LDA_DIVISOR_IN_BYTES ;
>>>> return
>>>> ( ld_ = bytesNbrRows / sizeof( _T ) ) ;
>>>> }
>>>>
>>>> size_t size1_ ;
>>>> size_t size2_ ;
>>>>
>>>> public:
>>>> typedef typename
>>>> _T value_type ;
>>>> typedef typename
>>>> ublas::matrix< _T, _L, _A > internal_matrix_type ;
>>>> typedef typename
>>>> ublas::matrix_range<
>>>> internal_matrix_type
>>>> > matrix_range_type ;
>>>> typedef typename
>>>> ublas::matrix_range<
>>>> const internal_matrix_type
>>>> > const_matrix_range_type ;
>>>>
>>>> typedef typename internal_matrix_type M ;
>>>> typedef typename M matrix_type;
>>>> typedef typename M::size_type size_type;
>>>> typedef typename M::difference_type difference_type;
>>>> typedef typename M::const_reference const_reference;
>>>> typedef typename boost::mpl::if_<boost::is_const<M>,
>>>> typename M::const_reference,
>>>> typename M::reference>::type reference;
>>>> typedef typename boost::mpl::if_<boost::is_const<M>,
>>>> typename M::const_closure_type,
>>>> typename M::closure_type>::type matrix_closure_type;
>>>> typedef typename ublas::basic_range<size_type, difference_type> 
>>>> range_type;
>>>>
>>>> typedef typename matrix_range_type self_type ;
>>>> typedef const self_type const_closure_type;
>>>> typedef self_type closure_type;
>>>> typedef typename ublas::storage_restrict_traits<typename 
>>>> M::storage_category,
>>>> ublas::dense_proxy_tag>::storage_category storage_category;
>>>> typedef typename M::orientation_category orientation_category;
>>>>
>>>> #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
>>>>           typedef indexed_iterator1<matrix_range<matrix_type>,
>>>>                                     typename 
>>>> subiterator1_type::iterator_category> iterator1;
>>>>           typedef indexed_iterator2<matrix_range<matrix_type>,
>>>>                                     typename 
>>>> subiterator2_type::iterator_category> iterator2;
>>>>           typedef indexed_const_iterator1<matrix_range<matrix_type>,
>>>>                                           typename 
>>>> const_subiterator1_type::iterator_category> const_iterator1;
>>>>           typedef indexed_const_iterator2<matrix_range<matrix_type>,
>>>>                                           typename 
>>>> const_subiterator2_type::iterator_category> const_iterator2;
>>>> #else
>>>>           class const_iterator1;
>>>>           class iterator1;
>>>>           class const_iterator2;
>>>>           class iterator2;
>>>> #endif
>>>>           typedef ublas::reverse_iterator_base1<const_iterator1> 
>>>> const_reverse_iterator1;
>>>>           typedef ublas::reverse_iterator_base1<iterator1> 
>>>> reverse_iterator1;
>>>>           typedef ublas::reverse_iterator_base2<const_iterator2> 
>>>> const_reverse_iterator2;
>>>>           typedef ublas::reverse_iterator_base2<iterator2> 
>>>> reverse_iterator2;
>>>> private :
>>>> typename internal_matrix_type matrix_ ;
>>>> public :
>>>> typedef MatrixGeneralT< _T, _L, _A> matrix_ge_self_type ;
>>>> MatrixGeneralT( bool optimalLDA = true )
>>>> :size1_(0), size2_(0)
>>>> {}
>>>> MatrixGeneralT( const size_t nbrRows,
>>>> const size_t nbrCols = 0,
>>>> bool optimalLDA = true
>>>> )
>>>> : size1_(nbrRows), size2_(nbrCols )
>>>> , matrix_( determineLDA( nbrRows, optimalLDA), nbrCols )
>>>> {}
>>>> virtual ~MatrixGeneralT()
>>>> {}
>>>>
>>>>
>>>> template <typename _me>
>>>> matrix_ge_self_type  &
>>>> operator  = (_me const & me )
>>>> {
>>>> matrix_range_type
>>>> (
>>>> matrix_
>>>> , ublas::range( 0, size1_ )
>>>> , ublas::range( 0, size2_ )
>>>> )  = me ;
>>>> return
>>>> *this ;
>>>> }
>>>> operator matrix_range_type()
>>>> {
>>>> return
>>>> matrix_range_type
>>>> (
>>>> matrix_
>>>> , ublas::range( 0, size1_ )
>>>> , ublas::range( 0, size2_ )
>>>> )
>>>> }
>>>> //operator const_matrix_range_type () const {
>>>> // return
>>>> // const_matrix_range_type
>>>> // (
>>>> // matrix_
>>>> // , ublas::range( 0, size1_ )
>>>> // , ublas::range( 0, size2_ )
>>>> // )
>>>> //}
>>>>
>>>> } ;
>>>>
>>>> int main() {
>>>>
>>>> typedef MatrixGeneralT< double  > matrix ;
>>>>
>>>> matrix m(3,3) ;
>>>> ublas::vector< double > x( 3 ) ;
>>>> x[0] = 1; x[1] = 2 ; x[2] = 3;
>>>> // matrix::matrix_range_type mr(m) ;
>>>> m = ( ublas::outer_prod( x,x ) ) ;
>>>>
>>>>
>>>> matrix m1( 3, 5 );
>>>> matrix m2( 5, 4 ) ;
>>>> matrix m3( 3, 4 ) ;
>>>>
>>>> m3 = ublas::prod( m1, m2 ) ;
>>>> return
>>>> 1;
>>>> }
>>>> _______________________________________________
>>>> ublas mailing list
>>>> ublas_at_[hidden]
>>>> http://listarchives.boost.org/mailman/listinfo.cgi/ublas
>>>> Sent to: athanasios.iliopoulos.ctr.gr_at_[hidden]
>>> _______________________________________________
>>> ublas mailing list
>>> ublas_at_[hidden]
>>> http://listarchives.boost.org/mailman/listinfo.cgi/ublas
>>> Sent to: pmamales_at_[hidden]
>> _______________________________________________
>> ublas mailing list
>> ublas_at_[hidden]
>> http://listarchives.boost.org/mailman/listinfo.cgi/ublas
>> Sent to: athanasios.iliopoulos.ctr.gr_at_[hidden]
>
> _______________________________________________
> ublas mailing list
> ublas_at_[hidden]
> http://listarchives.boost.org/mailman/listinfo.cgi/ublas
> Sent to: pmamales_at_[hidden]
>
> _______________________________________________
> ublas mailing list
> ublas_at_[hidden]
> http://listarchives.boost.org/mailman/listinfo.cgi/ublas
> Sent to: athanasios.iliopoulos.ctr.gr_at_[hidden]