|
|
Ublas : |
Subject: Re: [ublas] matrix_range for matrix
From: pmamales_at_[hidden]
Date: 2013-07-05 13:12:06
Naso,
Btw, just rechecked (debug) the outer_prod and works OK as is.
Petros
---- 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://lists.boost.org/mailman/listinfo.cgi/ublas
> > Sent to: athanasios.iliopoulos.ctr.gr_at_[hidden]
>
> _______________________________________________
> ublas mailing list
> ublas_at_[hidden]
> http://lists.boost.org/mailman/listinfo.cgi/ublas
> Sent to: pmamales_at_[hidden]