$include_dir="/home/hyper-archives/ublas/include"; include("$include_dir/msg-header.inc") ?>
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://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]