Subject: Re: [ublas] matrix_range for matrix
From: oswin krause (oswin.krause_at_[hidden])
Date: 2013-07-06 04:28:49
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: ublas_at_[hidden]
> 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: 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]
>
> _______________________________________________
> 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: Oswin.Krause_at_[hidden]