Re: matrix_range for matrix

Hi Oswin,
The problem with this construction, is that you cannot resize this matrix.
Petros

-----Original Message----- 
From: oswin krause 
Sent: Saturday, July 06, 2013 4:28 AM 
To: [email protected] 
Subject: Re: [ublas] matrix_range for matrix 

Hi,

there is no need to dig into the internals of ublas. In your case you 
could implment your class in terms of ublas::matrix and 
ublas::matrix_range and than just forward the public interface of the 
range. This could look similar to this (just the idea, no compilable 
code...):

class MatrixGeneralT{
public:
     MatrixGeneralT(
        const size_t nbrRows,
        const size_t nbrCols = 0
    )
        : matrix_( determineLDA( nbrRows ), nbrCols )
        , interface_(subrange(matrix_,0,nbrRows,0,nbrCols))//connect the 
matrix with the subrange
    {}
    //...
     typedef typename matrix_range_type::iterator1 iterator1;
     typedef typename matrix_range_type::iterator2 iterator2;
     //....
    const_reference operator()(size_type i, size_type j)const{
        return interface_(i,j);
    }
    reference operator()(size_type i, size_type j){
        return interface_(i,j);
    }
    size_type size1()const{
        return interface_.size1();
    }
    size_type size2()const{
        return interface_.size1();
    }

    //rest in the same manner: begin1/begin2 end1/end2 find1/find2 and 
that should be everything needed i guess...


private:
     internal_matrix_type matrix_;
     matrix_range_type interface_;
};

On 06.07.2013 00:34, Petros wrote:
> Naso,
>
> Other than a)the extra constructor with the default nbrCols argument 
> and the b)namespace alias vs.
> the keyword using ( actually the latter is much more loose and better 
> be avoided, since other things get
> defined as well ) the two pieces of code are identical.
> ( the c/tor uses size2_ before it is initialized - this is the problem 
> in icl as well, only icl does not
> make a fuss about it if you don't use it ;-) ).
>
> I am using the intel 12.1 compiler with msvc10 on windows and c++0x.
> I don't understand why gcc would have any issue (maybe forgot the 
> std0x flag?) w/it but I do
> not build with gcc, anyway..
>
> As far as rethinking the problem, I have thought about it a lot and 
> about a year ago as well and
> cannot see any other way, without getting into the ublas internals - 
> in the end of the day, I
> am only a user of ublas.
> The way I see it, I need my matrix for 2 things: its own public 
> interface - I can take care of it, and its
> connectivity with the ET of ublas. The latter is in the form of either 
> assignments ( again I can take
> care of it ) and how it behaves on the lhs of expressions. The prod 
> example, shows that the trick of the
> auto cast probably does the job.
>
> I know it is not ideal, since it is uncontrollable and could "collide" 
> with something else. I am hopeful because, in reality,
> it is no more than another matrix expression . I give it the look and 
> feel of a matrix_range but "package" together the
> matrix it refers to.
> Maybe somebody else, who will have a better idea, given the 
> constraints, will come up and present it
> in the forum, or even better add it to the code. Until then ..
>
> Thank you for all your help,
> Petros
>
>
> -----Original Message----- From: Nasos Iliopoulos
> Sent: Friday, July 05, 2013 4:44 PM
> To: [email protected]
> Subject: Re: [ublas] matrix_range for matrix
>
> 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: [email protected]
>> 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, [email protected] 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 <[email protected]> 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, [email protected] 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
>>>>> [email protected]
>>>>> http://lists.boost.org/mailman/listinfo.cgi/ublas
>>>>> Sent to: [email protected]
>>>> _______________________________________________
>>>> ublas mailing list
>>>> [email protected]
>>>> http://lists.boost.org/mailman/listinfo.cgi/ublas
>>>> Sent to: [email protected]
>>> _______________________________________________
>>> ublas mailing list
>>> [email protected]
>>> http://lists.boost.org/mailman/listinfo.cgi/ublas
>>> Sent to: [email protected]
>>
>> _______________________________________________
>> ublas mailing list
>> [email protected]
>> http://lists.boost.org/mailman/listinfo.cgi/ublas
>> Sent to: [email protected]
>>
>> _______________________________________________
>> ublas mailing list
>> [email protected]
>> http://lists.boost.org/mailman/listinfo.cgi/ublas
>> Sent to: [email protected]
>
> _______________________________________________
> ublas mailing list
> [email protected]
> http://lists.boost.org/mailman/listinfo.cgi/ublas
> Sent to: [email protected]
>
> _______________________________________________
> ublas mailing list
> [email protected]
> http://lists.boost.org/mailman/listinfo.cgi/ublas
> Sent to: [email protected]

_______________________________________________
ublas mailing list
[email protected]
http://lists.boost.org/mailman/listinfo.cgi/ublas
Sent to: [email protected]

_______________________________________________
ublas mailing list
[email protected]
http://lists.boost.org/mailman/listinfo.cgi/ublas
Sent to: [email protected]