SMatDMatMultExpr.h
Go to the documentation of this file.
1 //=================================================================================================
33 //=================================================================================================
34 
35 #ifndef _BLAZE_MATH_EXPRESSIONS_SMATDMATMULTEXPR_H_
36 #define _BLAZE_MATH_EXPRESSIONS_SMATDMATMULTEXPR_H_
37 
38 
39 //*************************************************************************************************
40 // Includes
41 //*************************************************************************************************
42 
43 #include <blaze/math/Aliases.h>
50 #include <blaze/math/Exception.h>
55 #include <blaze/math/Functions.h>
56 #include <blaze/math/shims/Reset.h>
58 #include <blaze/math/SIMD.h>
98 #include <blaze/util/Assert.h>
99 #include <blaze/util/EnableIf.h>
101 #include <blaze/util/InvalidType.h>
103 #include <blaze/util/mpl/And.h>
104 #include <blaze/util/mpl/If.h>
105 #include <blaze/util/mpl/Or.h>
106 #include <blaze/util/Types.h>
108 
109 
110 namespace blaze {
111 
112 //=================================================================================================
113 //
114 // CLASS SMATDMATMULTEXPR
115 //
116 //=================================================================================================
117 
118 //*************************************************************************************************
125 template< typename MT1 // Type of the left-hand side sparse matrix
126  , typename MT2 > // Type of the right-hand side dense matrix
127 class SMatDMatMultExpr : public DenseMatrix< SMatDMatMultExpr<MT1,MT2>, false >
128  , private MatMatMultExpr
129  , private Computation
130 {
131  private:
132  //**Type definitions****************************************************************************
139  //**********************************************************************************************
140 
141  //**********************************************************************************************
143  enum : bool { evaluateLeft = IsComputation<MT1>::value || RequiresEvaluation<MT1>::value };
144  //**********************************************************************************************
145 
146  //**********************************************************************************************
148  enum : bool { evaluateRight = IsComputation<MT2>::value || RequiresEvaluation<MT2>::value };
149  //**********************************************************************************************
150 
151  //**********************************************************************************************
153 
157  template< typename T1, typename T2, typename T3 >
158  struct IsEvaluationRequired {
159  enum : bool { value = ( evaluateLeft || evaluateRight ) };
160  };
162  //**********************************************************************************************
163 
164  //**********************************************************************************************
166 
169  template< typename T1, typename T2, typename T3 >
170  struct UseVectorizedKernel {
171  enum : bool { value = useOptimizedKernels &&
172  !IsDiagonal<T3>::value &&
173  T1::simdEnabled && T3::simdEnabled &&
174  IsRowMajorMatrix<T1>::value &&
175  AreSIMDCombinable< ElementType_<T1>
176  , ElementType_<T2>
177  , ElementType_<T3> >::value &&
178  HasSIMDAdd< ElementType_<T2>, ElementType_<T3> >::value &&
179  HasSIMDMult< ElementType_<T2>, ElementType_<T3> >::value };
180  };
182  //**********************************************************************************************
183 
184  //**********************************************************************************************
186 
190  template< typename T1, typename T2, typename T3 >
191  struct UseOptimizedKernel {
192  enum : bool { value = useOptimizedKernels &&
193  !UseVectorizedKernel<T1,T2,T3>::value &&
194  !IsDiagonal<T3>::value &&
195  !IsResizable< ElementType_<T1> >::value &&
196  !IsResizable<ET1>::value };
197  };
199  //**********************************************************************************************
200 
201  //**********************************************************************************************
203 
206  template< typename T1, typename T2, typename T3 >
207  struct UseDefaultKernel {
208  enum : bool { value = !UseVectorizedKernel<T1,T2,T3>::value &&
209  !UseOptimizedKernel<T1,T2,T3>::value };
210  };
212  //**********************************************************************************************
213 
214  public:
215  //**Type definitions****************************************************************************
222  typedef const ElementType ReturnType;
223  typedef const ResultType CompositeType;
224 
226  typedef If_< IsExpression<MT1>, const MT1, const MT1& > LeftOperand;
227 
229  typedef If_< IsExpression<MT2>, const MT2, const MT2& > RightOperand;
230 
233 
236  //**********************************************************************************************
237 
238  //**Compilation flags***************************************************************************
240  enum : bool { simdEnabled = !IsDiagonal<MT2>::value &&
241  MT2::simdEnabled &&
244 
246  enum : bool { smpAssignable = !evaluateLeft && MT1::smpAssignable &&
247  !evaluateRight && MT2::smpAssignable };
248  //**********************************************************************************************
249 
250  //**SIMD properties*****************************************************************************
252  enum : size_t { SIMDSIZE = SIMDTrait<ElementType>::size };
253  //**********************************************************************************************
254 
255  //**Constructor*********************************************************************************
261  explicit inline SMatDMatMultExpr( const MT1& lhs, const MT2& rhs ) noexcept
262  : lhs_( lhs ) // Left-hand side sparse matrix of the multiplication expression
263  , rhs_( rhs ) // Right-hand side dense matrix of the multiplication expression
264  {
265  BLAZE_INTERNAL_ASSERT( lhs.columns() == rhs.rows(), "Invalid matrix sizes" );
266  }
267  //**********************************************************************************************
268 
269  //**Access operator*****************************************************************************
276  inline ReturnType operator()( size_t i, size_t j ) const {
277  BLAZE_INTERNAL_ASSERT( i < lhs_.rows() , "Invalid row access index" );
278  BLAZE_INTERNAL_ASSERT( j < rhs_.columns(), "Invalid column access index" );
279 
280  if( IsDiagonal<MT1>::value ) {
281  return lhs_(i,i) * rhs_(i,j);
282  }
283  else if( IsDiagonal<MT2>::value ) {
284  return lhs_(i,j) * rhs_(j,j);
285  }
287  const size_t begin( ( IsUpper<MT1>::value )
288  ?( ( IsLower<MT2>::value )
289  ?( max( ( IsStrictlyUpper<MT1>::value ? i+1UL : i )
290  , ( IsStrictlyLower<MT2>::value ? j+1UL : j ) ) )
291  :( IsStrictlyUpper<MT1>::value ? i+1UL : i ) )
292  :( ( IsLower<MT2>::value )
293  ?( IsStrictlyLower<MT2>::value ? j+1UL : j )
294  :( 0UL ) ) );
295  const size_t end( ( IsLower<MT1>::value )
296  ?( ( IsUpper<MT2>::value )
297  ?( min( ( IsStrictlyLower<MT1>::value ? i : i+1UL )
298  , ( IsStrictlyUpper<MT2>::value ? j : j+1UL ) ) )
299  :( IsStrictlyLower<MT1>::value ? i : i+1UL ) )
300  :( ( IsUpper<MT2>::value )
301  ?( IsStrictlyUpper<MT2>::value ? j : j+1UL )
302  :( lhs_.columns() ) ) );
303 
304  if( begin >= end ) return ElementType();
305 
306  const size_t n( end - begin );
307 
308  return subvector( row( lhs_, i ), begin, n ) * subvector( column( rhs_, j ), begin, n );
309  }
310  else {
311  return row( lhs_, i ) * column( rhs_, j );
312  }
313  }
314  //**********************************************************************************************
315 
316  //**At function*********************************************************************************
324  inline ReturnType at( size_t i, size_t j ) const {
325  if( i >= lhs_.rows() ) {
326  BLAZE_THROW_OUT_OF_RANGE( "Invalid row access index" );
327  }
328  if( j >= rhs_.columns() ) {
329  BLAZE_THROW_OUT_OF_RANGE( "Invalid column access index" );
330  }
331  return (*this)(i,j);
332  }
333  //**********************************************************************************************
334 
335  //**Rows function*******************************************************************************
340  inline size_t rows() const noexcept {
341  return lhs_.rows();
342  }
343  //**********************************************************************************************
344 
345  //**Columns function****************************************************************************
350  inline size_t columns() const noexcept {
351  return rhs_.columns();
352  }
353  //**********************************************************************************************
354 
355  //**Left operand access*************************************************************************
360  inline LeftOperand leftOperand() const noexcept {
361  return lhs_;
362  }
363  //**********************************************************************************************
364 
365  //**Right operand access************************************************************************
370  inline RightOperand rightOperand() const noexcept {
371  return rhs_;
372  }
373  //**********************************************************************************************
374 
375  //**********************************************************************************************
381  template< typename T >
382  inline bool canAlias( const T* alias ) const noexcept {
383  return ( lhs_.isAliased( alias ) || rhs_.isAliased( alias ) );
384  }
385  //**********************************************************************************************
386 
387  //**********************************************************************************************
393  template< typename T >
394  inline bool isAliased( const T* alias ) const noexcept {
395  return ( lhs_.isAliased( alias ) || rhs_.isAliased( alias ) );
396  }
397  //**********************************************************************************************
398 
399  //**********************************************************************************************
404  inline bool isAligned() const noexcept {
405  return rhs_.isAligned();
406  }
407  //**********************************************************************************************
408 
409  //**********************************************************************************************
414  inline bool canSMPAssign() const noexcept {
415  return ( rows() * columns() >= SMP_SMATDMATMULT_THRESHOLD );
416  }
417  //**********************************************************************************************
418 
419  private:
420  //**Member variables****************************************************************************
421  LeftOperand lhs_;
422  RightOperand rhs_;
423  //**********************************************************************************************
424 
425  //**Assignment to dense matrices****************************************************************
438  template< typename MT // Type of the target dense matrix
439  , bool SO > // Storage order of the target dense matrix
440  friend inline void assign( DenseMatrix<MT,SO>& lhs, const SMatDMatMultExpr& rhs )
441  {
443 
444  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
445  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
446 
447  LT A( serial( rhs.lhs_ ) ); // Evaluation of the left-hand side sparse matrix operand
448  RT B( serial( rhs.rhs_ ) ); // Evaluation of the right-hand side dense matrix operand
449 
450  BLAZE_INTERNAL_ASSERT( A.rows() == rhs.lhs_.rows() , "Invalid number of rows" );
451  BLAZE_INTERNAL_ASSERT( A.columns() == rhs.lhs_.columns(), "Invalid number of columns" );
452  BLAZE_INTERNAL_ASSERT( B.rows() == rhs.rhs_.rows() , "Invalid number of rows" );
453  BLAZE_INTERNAL_ASSERT( B.columns() == rhs.rhs_.columns(), "Invalid number of columns" );
454  BLAZE_INTERNAL_ASSERT( A.rows() == (~lhs).rows() , "Invalid number of rows" );
455  BLAZE_INTERNAL_ASSERT( B.columns() == (~lhs).columns() , "Invalid number of columns" );
456 
457  SMatDMatMultExpr::selectAssignKernel( ~lhs, A, B );
458  }
460  //**********************************************************************************************
461 
462  //**Default assignment to dense matrices********************************************************
476  template< typename MT3 // Type of the left-hand side target matrix
477  , typename MT4 // Type of the left-hand side matrix operand
478  , typename MT5 > // Type of the right-hand side matrix operand
480  selectAssignKernel( MT3& C, const MT4& A, const MT5& B )
481  {
483 
484  const size_t block( Or< IsRowMajorMatrix<MT3>, IsDiagonal<MT5> >::value ? B.columns() : 64UL );
485 
486  reset( C );
487 
488  for( size_t jj=0UL; jj<B.columns(); jj+=block )
489  {
490  const size_t jtmp( min( jj+block, B.columns() ) );
491 
492  for( size_t i=0UL; i<A.rows(); ++i )
493  {
494  ConstIterator element( A.begin(i) );
495  const ConstIterator end( A.end(i) );
496 
497  for( ; element!=end; ++element )
498  {
499  const size_t i1( element->index() );
500 
502  {
503  C(i,i1) = element->value() * B(i1,i1);
504  }
505  else
506  {
507  const size_t jbegin( ( IsUpper<MT5>::value )
508  ?( max( jj, ( IsStrictlyUpper<MT5>::value ? i1+1UL : i1 ) ) )
509  :( jj ) );
510  const size_t jend( ( IsLower<MT5>::value )
511  ?( min( jtmp, ( IsStrictlyLower<MT5>::value ? i1 : i1+1UL ) ) )
512  :( jtmp ) );
513 
514  if( IsTriangular<MT5>::value && jbegin >= jend )
515  continue;
516 
517  BLAZE_INTERNAL_ASSERT( jbegin <= jend, "Invalid loop indices detected" );
518 
519  for( size_t j=jbegin; j<jend; ++j ) {
520  if( isDefault( C(i,j) ) )
521  C(i,j) = element->value() * B(i1,j);
522  else
523  C(i,j) += element->value() * B(i1,j);
524  }
525  }
526  }
527  }
528  }
529  }
531  //**********************************************************************************************
532 
533  //**Optimized assignment to dense matrices******************************************************
547  template< typename MT3 // Type of the left-hand side target matrix
548  , typename MT4 // Type of the left-hand side matrix operand
549  , typename MT5 > // Type of the right-hand side matrix operand
550  static inline EnableIf_< UseOptimizedKernel<MT3,MT4,MT5> >
551  selectAssignKernel( MT3& C, const MT4& A, const MT5& B )
552  {
553  typedef ConstIterator_<MT4> ConstIterator;
554 
555  const size_t block( IsRowMajorMatrix<MT3>::value ? B.columns() : 64UL );
556 
557  reset( C );
558 
559  for( size_t jj=0UL; jj<B.columns(); jj+=block )
560  {
561  const size_t jtmp( min( jj+block, B.columns() ) );
562 
563  for( size_t i=0UL; i<A.rows(); ++i )
564  {
565  const ConstIterator end( A.end(i) );
566  ConstIterator element( A.begin(i) );
567 
568  const size_t nonzeros( A.nonZeros(i) );
569  const size_t kpos( nonzeros & size_t(-4) );
570  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
571 
572  for( size_t k=0UL; k<kpos; k+=4UL )
573  {
574  const size_t i1( element->index() );
575  const ET1 v1( element->value() );
576  ++element;
577  const size_t i2( element->index() );
578  const ET1 v2( element->value() );
579  ++element;
580  const size_t i3( element->index() );
581  const ET1 v3( element->value() );
582  ++element;
583  const size_t i4( element->index() );
584  const ET1 v4( element->value() );
585  ++element;
586 
587  BLAZE_INTERNAL_ASSERT( i1 < i2 && i2 < i3 && i3 < i4, "Invalid sparse matrix index detected" );
588 
589  const size_t jbegin( ( IsUpper<MT5>::value )
590  ?( max( jj, ( IsStrictlyUpper<MT5>::value ? i1+1UL : i1 ) ) )
591  :( jj ) );
592  const size_t jend( ( IsLower<MT5>::value )
593  ?( min( jtmp, ( IsStrictlyLower<MT5>::value ? i4 : i4+1UL ) ) )
594  :( jtmp ) );
595 
596  if( IsTriangular<MT5>::value && jbegin >= jend )
597  continue;
598 
599  BLAZE_INTERNAL_ASSERT( jbegin <= jend, "Invalid loop indices detected" );
600 
601  const size_t jnum( jend - jbegin );
602  const size_t jpos( jbegin + ( jnum & size_t(-4) ) );
603  BLAZE_INTERNAL_ASSERT( ( jbegin + jnum - ( jnum % 4UL ) ) == jpos, "Invalid end calculation" );
604 
605  for( size_t j=jbegin; j<jpos; j+=4UL ) {
606  C(i,j ) += v1 * B(i1,j ) + v2 * B(i2,j ) + v3 * B(i3,j ) + v4 * B(i4,j );
607  C(i,j+1UL) += v1 * B(i1,j+1UL) + v2 * B(i2,j+1UL) + v3 * B(i3,j+1UL) + v4 * B(i4,j+1UL);
608  C(i,j+2UL) += v1 * B(i1,j+2UL) + v2 * B(i2,j+2UL) + v3 * B(i3,j+2UL) + v4 * B(i4,j+2UL);
609  C(i,j+3UL) += v1 * B(i1,j+3UL) + v2 * B(i2,j+3UL) + v3 * B(i3,j+3UL) + v4 * B(i4,j+3UL);
610  }
611  for( size_t j=jpos; j<jend; ++j ) {
612  C(i,j) += v1 * B(i1,j) + v2 * B(i2,j) + v3 * B(i3,j) + v4 * B(i4,j);
613  }
614  }
615 
616  for( ; element!=end; ++element )
617  {
618  const size_t i1( element->index() );
619  const ET1 v1( element->value() );
620 
621  const size_t jbegin( ( IsUpper<MT5>::value )
622  ?( max( jj, ( IsStrictlyUpper<MT5>::value ? i1+1UL : i1 ) ) )
623  :( jj ) );
624  const size_t jend( ( IsLower<MT5>::value )
625  ?( min( jtmp, ( IsStrictlyLower<MT5>::value ? i1 : i1+1UL ) ) )
626  :( jtmp ) );
627 
628  if( IsTriangular<MT5>::value && jbegin >= jend )
629  continue;
630 
631  BLAZE_INTERNAL_ASSERT( jbegin <= jend, "Invalid loop indices detected" );
632 
633  const size_t jnum( jend - jbegin );
634  const size_t jpos( jbegin + ( jnum & size_t(-4) ) );
635  BLAZE_INTERNAL_ASSERT( ( jbegin + jnum - ( jnum % 4UL ) ) == jpos, "Invalid end calculation" );
636 
637  for( size_t j=jbegin; j<jpos; j+=4UL ) {
638  C(i,j ) += v1 * B(i1,j );
639  C(i,j+1UL) += v1 * B(i1,j+1UL);
640  C(i,j+2UL) += v1 * B(i1,j+2UL);
641  C(i,j+3UL) += v1 * B(i1,j+3UL);
642  }
643  for( size_t j=jpos; j<jend; ++j ) {
644  C(i,j) += v1 * B(i1,j);
645  }
646  }
647  }
648  }
649  }
651  //**********************************************************************************************
652 
653  //**Vectorized assignment to dense matrices*****************************************************
667  template< typename MT3 // Type of the left-hand side target matrix
668  , typename MT4 // Type of the left-hand side matrix operand
669  , typename MT5 > // Type of the right-hand side matrix operand
670  static inline EnableIf_< UseVectorizedKernel<MT3,MT4,MT5> >
671  selectAssignKernel( MT3& C, const MT4& A, const MT5& B )
672  {
673  typedef ConstIterator_<MT4> ConstIterator;
674 
675  const bool remainder( !IsPadded<MT3>::value || !IsPadded<MT5>::value );
676 
677  reset( C );
678 
679  for( size_t i=0UL; i<A.rows(); ++i )
680  {
681  const ConstIterator end( A.end(i) );
682  ConstIterator element( A.begin(i) );
683 
684  const size_t nonzeros( A.nonZeros(i) );
685  const size_t kpos( nonzeros & size_t(-4) );
686  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
687 
688  for( size_t k=0UL; k<kpos; k+=4UL )
689  {
690  const size_t i1( element->index() );
691  const ET1 v1( element->value() );
692  ++element;
693  const size_t i2( element->index() );
694  const ET1 v2( element->value() );
695  ++element;
696  const size_t i3( element->index() );
697  const ET1 v3( element->value() );
698  ++element;
699  const size_t i4( element->index() );
700  const ET1 v4( element->value() );
701  ++element;
702 
703  BLAZE_INTERNAL_ASSERT( i1 < i2 && i2 < i3 && i3 < i4, "Invalid sparse matrix index detected" );
704 
705  const SIMDType xmm1( set( v1 ) );
706  const SIMDType xmm2( set( v2 ) );
707  const SIMDType xmm3( set( v3 ) );
708  const SIMDType xmm4( set( v4 ) );
709 
710  const size_t jbegin( ( IsUpper<MT5>::value )
711  ?( ( IsStrictlyUpper<MT5>::value ? i1+1UL : i1 ) & size_t(-SIMDSIZE) )
712  :( 0UL ) );
713  const size_t jend( ( IsLower<MT5>::value )
714  ?( IsStrictlyLower<MT5>::value ? i4 : i4+1UL )
715  :( B.columns() ) );
716  BLAZE_INTERNAL_ASSERT( jbegin <= jend, "Invalid loop indices detected" );
717 
718  const size_t jpos( remainder ? ( jend & size_t(-SIMDSIZE) ) : jend );
719  BLAZE_INTERNAL_ASSERT( !remainder || ( jend - ( jend % (SIMDSIZE) ) ) == jpos, "Invalid end calculation" );
720 
721  size_t j( jbegin );
722 
723  for( ; j<jpos; j+=SIMDSIZE ) {
724  C.store( i, j, C.load(i,j) + xmm1 * B.load(i1,j) + xmm2 * B.load(i2,j) + xmm3 * B.load(i3,j) + xmm4 * B.load(i4,j) );
725  }
726  for( ; remainder && j<jend; ++j ) {
727  C(i,j) += v1 * B(i1,j) + v2 * B(i2,j) + v3 * B(i3,j) + v4 * B(i4,j);
728  }
729  }
730 
731  for( ; element!=end; ++element )
732  {
733  const size_t i1( element->index() );
734  const ET1 v1( element->value() );
735 
736  const SIMDType xmm1( set( v1 ) );
737 
738  const size_t jbegin( ( IsUpper<MT5>::value )
739  ?( ( IsStrictlyUpper<MT5>::value ? i1+1UL : i1 ) & size_t(-SIMDSIZE) )
740  :( 0UL ) );
741  const size_t jend( ( IsLower<MT5>::value )
742  ?( IsStrictlyLower<MT5>::value ? i1 : i1+1UL )
743  :( B.columns() ) );
744  BLAZE_INTERNAL_ASSERT( jbegin <= jend, "Invalid loop indices detected" );
745 
746  const size_t jpos( remainder ? ( jend & size_t(-SIMDSIZE) ) : jend );
747  BLAZE_INTERNAL_ASSERT( !remainder || ( jend - ( jend % (SIMDSIZE) ) ) == jpos, "Invalid end calculation" );
748 
749  size_t j( jbegin );
750 
751  for( ; j<jpos; j+=SIMDSIZE ) {
752  C.store( i, j, C.load(i,j) + xmm1 * B.load(i1,j) );
753  }
754  for( ; remainder && j<jend; ++j ) {
755  C(i,j) += v1 * B(i1,j);
756  }
757  }
758  }
759  }
761  //**********************************************************************************************
762 
763  //**Assignment to sparse matrices***************************************************************
776  template< typename MT // Type of the target sparse matrix
777  , bool SO > // Storage order of the target sparse matrix
778  friend inline void assign( SparseMatrix<MT,SO>& lhs, const SMatDMatMultExpr& rhs )
779  {
781 
782  typedef IfTrue_< SO, OppositeType, ResultType > TmpType;
783 
789  BLAZE_CONSTRAINT_MUST_BE_REFERENCE_TYPE( CompositeType_<TmpType> );
790 
791  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
792  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
793 
794  const TmpType tmp( serial( rhs ) );
795  assign( ~lhs, tmp );
796  }
798  //**********************************************************************************************
799 
800  //**Addition assignment to dense matrices*******************************************************
813  template< typename MT // Type of the target dense matrix
814  , bool SO > // Storage order of the target dense matrix
815  friend inline void addAssign( DenseMatrix<MT,SO>& lhs, const SMatDMatMultExpr& rhs )
816  {
818 
819  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
820  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
821 
822  LT A( serial( rhs.lhs_ ) ); // Evaluation of the left-hand side sparse matrix operand
823  RT B( serial( rhs.rhs_ ) ); // Evaluation of the right-hand side dense matrix operand
824 
825  BLAZE_INTERNAL_ASSERT( A.rows() == rhs.lhs_.rows() , "Invalid number of rows" );
826  BLAZE_INTERNAL_ASSERT( A.columns() == rhs.lhs_.columns(), "Invalid number of columns" );
827  BLAZE_INTERNAL_ASSERT( B.rows() == rhs.rhs_.rows() , "Invalid number of rows" );
828  BLAZE_INTERNAL_ASSERT( B.columns() == rhs.rhs_.columns(), "Invalid number of columns" );
829  BLAZE_INTERNAL_ASSERT( A.rows() == (~lhs).rows() , "Invalid number of rows" );
830  BLAZE_INTERNAL_ASSERT( B.columns() == (~lhs).columns() , "Invalid number of columns" );
831 
832  SMatDMatMultExpr::selectAddAssignKernel( ~lhs, A, B );
833  }
835  //**********************************************************************************************
836 
837  //**Default addition assignment to dense matrices***********************************************
851  template< typename MT3 // Type of the left-hand side target matrix
852  , typename MT4 // Type of the left-hand side matrix operand
853  , typename MT5 > // Type of the right-hand side matrix operand
854  static inline EnableIf_< UseDefaultKernel<MT3,MT4,MT5> >
855  selectAddAssignKernel( MT3& C, const MT4& A, const MT5& B )
856  {
857  typedef ConstIterator_<MT4> ConstIterator;
858 
859  const size_t block( Or< IsRowMajorMatrix<MT3>, IsDiagonal<MT5> >::value ? B.columns() : 64UL );
860 
861  for( size_t jj=0UL; jj<B.columns(); jj+=block )
862  {
863  const size_t jtmp( min( jj+block, B.columns() ) );
864 
865  for( size_t i=0UL; i<A.rows(); ++i )
866  {
867  const ConstIterator end( A.end(i) );
868  ConstIterator element( A.begin(i) );
869 
870  for( ; element!=end; ++element )
871  {
872  const size_t i1( element->index() );
873 
874  if( IsDiagonal<MT5>::value )
875  {
876  C(i,i1) += element->value() * B(i1,i1);
877  }
878  else
879  {
880  const size_t jbegin( ( IsUpper<MT5>::value )
881  ?( max( jj, ( IsStrictlyUpper<MT5>::value ? i1+1UL : i1 ) ) )
882  :( jj ) );
883  const size_t jend( ( IsLower<MT5>::value )
884  ?( min( jtmp, ( IsStrictlyLower<MT5>::value ? i1 : i1+1UL ) ) )
885  :( jtmp ) );
886 
887  if( IsTriangular<MT5>::value && jbegin >= jend )
888  continue;
889 
890  BLAZE_INTERNAL_ASSERT( jbegin <= jend, "Invalid loop indices detected" );
891 
892  const size_t jnum( jend - jbegin );
893  const size_t jpos( jbegin + ( jnum & size_t(-4) ) );
894  BLAZE_INTERNAL_ASSERT( ( jbegin + jnum - ( jnum % 4UL ) ) == jpos, "Invalid end calculation" );
895 
896  for( size_t j=jbegin; j<jpos; j+=4UL ) {
897  C(i,j ) += element->value() * B(i1,j );
898  C(i,j+1UL) += element->value() * B(i1,j+1UL);
899  C(i,j+2UL) += element->value() * B(i1,j+2UL);
900  C(i,j+3UL) += element->value() * B(i1,j+3UL);
901  }
902  for( size_t j=jpos; j<jend; ++j ) {
903  C(i,j) += element->value() * B(i1,j);
904  }
905  }
906  }
907  }
908  }
909  }
911  //**********************************************************************************************
912 
913  //**Optimized addition assignment to dense matrices*********************************************
927  template< typename MT3 // Type of the left-hand side target matrix
928  , typename MT4 // Type of the left-hand side matrix operand
929  , typename MT5 > // Type of the right-hand side matrix operand
930  static inline EnableIf_< UseOptimizedKernel<MT3,MT4,MT5> >
931  selectAddAssignKernel( MT3& C, const MT4& A, const MT5& B )
932  {
933  typedef ConstIterator_<MT4> ConstIterator;
934 
935  const size_t block( IsRowMajorMatrix<MT3>::value ? B.columns() : 64UL );
936 
937  for( size_t jj=0UL; jj<B.columns(); jj+=block )
938  {
939  const size_t jtmp( min( jj+block, B.columns() ) );
940 
941  for( size_t i=0UL; i<A.rows(); ++i )
942  {
943  const ConstIterator end( A.end(i) );
944  ConstIterator element( A.begin(i) );
945 
946  const size_t nonzeros( A.nonZeros(i) );
947  const size_t kpos( nonzeros & size_t(-4) );
948  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
949 
950  for( size_t k=0UL; k<kpos; k+=4UL )
951  {
952  const size_t i1( element->index() );
953  const ET1 v1( element->value() );
954  ++element;
955  const size_t i2( element->index() );
956  const ET1 v2( element->value() );
957  ++element;
958  const size_t i3( element->index() );
959  const ET1 v3( element->value() );
960  ++element;
961  const size_t i4( element->index() );
962  const ET1 v4( element->value() );
963  ++element;
964 
965  BLAZE_INTERNAL_ASSERT( i1 < i2 && i2 < i3 && i3 < i4, "Invalid sparse matrix index detected" );
966 
967  const size_t jbegin( ( IsUpper<MT5>::value )
968  ?( max( jj, ( IsStrictlyUpper<MT5>::value ? i1+1UL : i1 ) ) )
969  :( jj ) );
970  const size_t jend( ( IsLower<MT5>::value )
971  ?( min( jtmp, ( IsStrictlyLower<MT5>::value ? i4 : i4+1UL ) ) )
972  :( jtmp ) );
973 
974  if( IsTriangular<MT5>::value && jbegin >= jend )
975  continue;
976 
977  BLAZE_INTERNAL_ASSERT( jbegin <= jend, "Invalid loop indices detected" );
978 
979  const size_t jnum( jend - jbegin );
980  const size_t jpos( jbegin + ( jnum & size_t(-4) ) );
981  BLAZE_INTERNAL_ASSERT( ( jbegin + jnum - ( jnum % 4UL ) ) == jpos, "Invalid end calculation" );
982 
983  for( size_t j=jbegin; j<jpos; j+=4UL ) {
984  C(i,j ) += v1 * B(i1,j ) + v2 * B(i2,j ) + v3 * B(i3,j ) + v4 * B(i4,j );
985  C(i,j+1UL) += v1 * B(i1,j+1UL) + v2 * B(i2,j+1UL) + v3 * B(i3,j+1UL) + v4 * B(i4,j+1UL);
986  C(i,j+2UL) += v1 * B(i1,j+2UL) + v2 * B(i2,j+2UL) + v3 * B(i3,j+2UL) + v4 * B(i4,j+2UL);
987  C(i,j+3UL) += v1 * B(i1,j+3UL) + v2 * B(i2,j+3UL) + v3 * B(i3,j+3UL) + v4 * B(i4,j+3UL);
988  }
989  for( size_t j=jpos; j<jend; ++j ) {
990  C(i,j) += v1 * B(i1,j) + v2 * B(i2,j) + v3 * B(i3,j) + v4 * B(i4,j);
991  }
992  }
993 
994  for( ; element!=end; ++element )
995  {
996  const size_t i1( element->index() );
997  const ET1 v1( element->value() );
998 
999  const size_t jbegin( ( IsUpper<MT5>::value )
1000  ?( max( jj, ( IsStrictlyUpper<MT5>::value ? i1+1UL : i1 ) ) )
1001  :( jj ) );
1002  const size_t jend( ( IsLower<MT5>::value )
1003  ?( min( jtmp, ( IsStrictlyLower<MT5>::value ? i1 : i1+1UL ) ) )
1004  :( jtmp ) );
1005 
1006  if( IsTriangular<MT5>::value && jbegin >= jend )
1007  continue;
1008 
1009  BLAZE_INTERNAL_ASSERT( jbegin <= jend, "Invalid loop indices detected" );
1010 
1011  const size_t jnum( jend - jbegin );
1012  const size_t jpos( jbegin + ( jnum & size_t(-4) ) );
1013  BLAZE_INTERNAL_ASSERT( ( jbegin + jnum - ( jnum % 4UL ) ) == jpos, "Invalid end calculation" );
1014 
1015  for( size_t j=jbegin; j<jpos; j+=4UL ) {
1016  C(i,j ) += v1 * B(i1,j );
1017  C(i,j+1UL) += v1 * B(i1,j+1UL);
1018  C(i,j+2UL) += v1 * B(i1,j+2UL);
1019  C(i,j+3UL) += v1 * B(i1,j+3UL);
1020  }
1021  for( size_t j=jpos; j<jend; ++j ) {
1022  C(i,j) += v1 * B(i1,j);
1023  }
1024  }
1025  }
1026  }
1027  }
1029  //**********************************************************************************************
1030 
1031  //**Vectorized addition assignment to dense matrices********************************************
1045  template< typename MT3 // Type of the left-hand side target matrix
1046  , typename MT4 // Type of the left-hand side matrix operand
1047  , typename MT5 > // Type of the right-hand side matrix operand
1048  static inline EnableIf_< UseVectorizedKernel<MT3,MT4,MT5> >
1049  selectAddAssignKernel( MT3& C, const MT4& A, const MT5& B )
1050  {
1051  typedef ConstIterator_<MT4> ConstIterator;
1052 
1053  const bool remainder( !IsPadded<MT3>::value || !IsPadded<MT5>::value );
1054 
1055  for( size_t i=0UL; i<A.rows(); ++i )
1056  {
1057  const ConstIterator end( A.end(i) );
1058  ConstIterator element( A.begin(i) );
1059 
1060  const size_t nonzeros( A.nonZeros(i) );
1061  const size_t kpos( nonzeros & size_t(-4) );
1062  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
1063 
1064  for( size_t k=0UL; k<kpos; k+=4UL )
1065  {
1066  const size_t i1( element->index() );
1067  const ET1 v1( element->value() );
1068  ++element;
1069  const size_t i2( element->index() );
1070  const ET1 v2( element->value() );
1071  ++element;
1072  const size_t i3( element->index() );
1073  const ET1 v3( element->value() );
1074  ++element;
1075  const size_t i4( element->index() );
1076  const ET1 v4( element->value() );
1077  ++element;
1078 
1079  BLAZE_INTERNAL_ASSERT( i1 < i2 && i2 < i3 && i3 < i4, "Invalid sparse matrix index detected" );
1080 
1081  const SIMDType xmm1( set( v1 ) );
1082  const SIMDType xmm2( set( v2 ) );
1083  const SIMDType xmm3( set( v3 ) );
1084  const SIMDType xmm4( set( v4 ) );
1085 
1086  const size_t jbegin( ( IsUpper<MT5>::value )
1087  ?( ( IsStrictlyUpper<MT5>::value ? i1+1UL : i1 ) & size_t(-SIMDSIZE) )
1088  :( 0UL ) );
1089  const size_t jend( ( IsLower<MT5>::value )
1090  ?( IsStrictlyLower<MT5>::value ? i4 : i4+1UL )
1091  :( B.columns() ) );
1092  BLAZE_INTERNAL_ASSERT( jbegin <= jend, "Invalid loop indices detected" );
1093 
1094  const size_t jpos( remainder ? ( jend & size_t(-SIMDSIZE) ) : jend );
1095  BLAZE_INTERNAL_ASSERT( !remainder || ( jend - ( jend % (SIMDSIZE) ) ) == jpos, "Invalid end calculation" );
1096 
1097  size_t j( jbegin );
1098 
1099  for( ; j<jpos; j+=SIMDSIZE ) {
1100  C.store( i, j, C.load(i,j) + xmm1 * B.load(i1,j) + xmm2 * B.load(i2,j) + xmm3 * B.load(i3,j) + xmm4 * B.load(i4,j) );
1101  }
1102  for( ; remainder && j<jend; ++j ) {
1103  C(i,j) += v1 * B(i1,j) + v2 * B(i2,j) + v3 * B(i3,j) + v4 * B(i4,j);
1104  }
1105  }
1106 
1107  for( ; element!=end; ++element )
1108  {
1109  const size_t i1( element->index() );
1110  const ET1 v1( element->value() );
1111 
1112  const SIMDType xmm1( set( v1 ) );
1113 
1114  const size_t jbegin( ( IsUpper<MT5>::value )
1115  ?( ( IsStrictlyUpper<MT5>::value ? i1+1UL : i1 ) & size_t(-SIMDSIZE) )
1116  :( 0UL ) );
1117  const size_t jend( ( IsLower<MT5>::value )
1118  ?( IsStrictlyLower<MT5>::value ? i1 : i1+1UL )
1119  :( B.columns() ) );
1120  BLAZE_INTERNAL_ASSERT( jbegin <= jend, "Invalid loop indices detected" );
1121 
1122  const size_t jpos( remainder ? ( jend & size_t(-SIMDSIZE) ) : jend );
1123  BLAZE_INTERNAL_ASSERT( !remainder || ( jend - ( jend % (SIMDSIZE) ) ) == jpos, "Invalid end calculation" );
1124 
1125  size_t j( jbegin );
1126 
1127  for( ; j<jpos; j+=SIMDSIZE ) {
1128  C.store( i, j, C.load(i,j) + xmm1 * B.load(i1,j) );
1129  }
1130  for( ; remainder && j<jend; ++j ) {
1131  C(i,j) += v1 * B(i1,j);
1132  }
1133  }
1134  }
1135  }
1137  //**********************************************************************************************
1138 
1139  //**Addition assignment to sparse matrices******************************************************
1140  // No special implementation for the addition assignment to sparse matrices.
1141  //**********************************************************************************************
1142 
1143  //**Subtraction assignment to dense matrices****************************************************
1156  template< typename MT // Type of the target dense matrix
1157  , bool SO > // Storage order of the target dense matrix
1158  friend inline void subAssign( DenseMatrix<MT,SO>& lhs, const SMatDMatMultExpr& rhs )
1159  {
1161 
1162  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1163  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1164 
1165  LT A( serial( rhs.lhs_ ) ); // Evaluation of the left-hand side sparse matrix operand
1166  RT B( serial( rhs.rhs_ ) ); // Evaluation of the right-hand side dense matrix operand
1167 
1168  BLAZE_INTERNAL_ASSERT( A.rows() == rhs.lhs_.rows() , "Invalid number of rows" );
1169  BLAZE_INTERNAL_ASSERT( A.columns() == rhs.lhs_.columns(), "Invalid number of columns" );
1170  BLAZE_INTERNAL_ASSERT( B.rows() == rhs.rhs_.rows() , "Invalid number of rows" );
1171  BLAZE_INTERNAL_ASSERT( B.columns() == rhs.rhs_.columns(), "Invalid number of columns" );
1172  BLAZE_INTERNAL_ASSERT( A.rows() == (~lhs).rows() , "Invalid number of rows" );
1173  BLAZE_INTERNAL_ASSERT( B.columns() == (~lhs).columns() , "Invalid number of columns" );
1174 
1175  SMatDMatMultExpr::selectSubAssignKernel( ~lhs, A, B );
1176  }
1178  //**********************************************************************************************
1179 
1180  //**Default subtraction assignment to dense matrices********************************************
1194  template< typename MT3 // Type of the left-hand side target matrix
1195  , typename MT4 // Type of the left-hand side matrix operand
1196  , typename MT5 > // Type of the right-hand side matrix operand
1197  static inline EnableIf_< UseDefaultKernel<MT3,MT4,MT5> >
1198  selectSubAssignKernel( MT3& C, const MT4& A, const MT5& B )
1199  {
1200  typedef ConstIterator_<MT4> ConstIterator;
1201 
1202  const size_t block( Or< IsRowMajorMatrix<MT3>, IsDiagonal<MT5> >::value ? B.columns() : 64UL );
1203 
1204  for( size_t jj=0UL; jj<B.columns(); jj+=block )
1205  {
1206  const size_t jtmp( min( jj+block, B.columns() ) );
1207 
1208  for( size_t i=0UL; i<A.rows(); ++i )
1209  {
1210  const ConstIterator end( A.end(i) );
1211  ConstIterator element( A.begin(i) );
1212 
1213  for( ; element!=end; ++element )
1214  {
1215  const size_t i1( element->index() );
1216 
1217  if( IsDiagonal<MT5>::value )
1218  {
1219  C(i,i1) -= element->value() * B(i1,i1);
1220  }
1221  else
1222  {
1223  const size_t jbegin( ( IsUpper<MT5>::value )
1224  ?( max( jj, ( IsStrictlyUpper<MT5>::value ? i1+1UL : i1 ) ) )
1225  :( jj ) );
1226  const size_t jend( ( IsLower<MT5>::value )
1227  ?( min( jtmp, ( IsStrictlyLower<MT5>::value ? i1 : i1+1UL ) ) )
1228  :( jtmp ) );
1229 
1230  if( IsTriangular<MT5>::value && jbegin >= jend )
1231  continue;
1232 
1233  BLAZE_INTERNAL_ASSERT( jbegin <= jend, "Invalid loop indices detected" );
1234 
1235  const size_t jnum( jend - jbegin );
1236  const size_t jpos( jbegin + ( jnum & size_t(-4) ) );
1237  BLAZE_INTERNAL_ASSERT( ( jbegin + jnum - ( jnum % 4UL ) ) == jpos, "Invalid end calculation" );
1238 
1239  for( size_t j=jbegin; j<jpos; j+=4UL ) {
1240  C(i,j ) -= element->value() * B(i1,j );
1241  C(i,j+1UL) -= element->value() * B(i1,j+1UL);
1242  C(i,j+2UL) -= element->value() * B(i1,j+2UL);
1243  C(i,j+3UL) -= element->value() * B(i1,j+3UL);
1244  }
1245  for( size_t j=jpos; j<jend; ++j ) {
1246  C(i,j) -= element->value() * B(i1,j);
1247  }
1248  }
1249  }
1250  }
1251  }
1252  }
1254  //**********************************************************************************************
1255 
1256  //**Optimized subtraction assignment to dense matrices******************************************
1270  template< typename MT3 // Type of the left-hand side target matrix
1271  , typename MT4 // Type of the left-hand side matrix operand
1272  , typename MT5 > // Type of the right-hand side matrix operand
1273  static inline EnableIf_< UseOptimizedKernel<MT3,MT4,MT5> >
1274  selectSubAssignKernel( MT3& C, const MT4& A, const MT5& B )
1275  {
1276  typedef ConstIterator_<MT4> ConstIterator;
1277 
1278  const size_t block( IsRowMajorMatrix<MT3>::value ? B.columns() : 64UL );
1279 
1280  for( size_t jj=0UL; jj<B.columns(); jj+=block )
1281  {
1282  const size_t jtmp( min( jj+block, B.columns() ) );
1283 
1284  for( size_t i=0UL; i<A.rows(); ++i )
1285  {
1286  const ConstIterator end( A.end(i) );
1287  ConstIterator element( A.begin(i) );
1288 
1289  const size_t nonzeros( A.nonZeros(i) );
1290  const size_t kpos( nonzeros & size_t(-4) );
1291  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
1292 
1293  for( size_t k=0UL; k<kpos; k+=4UL )
1294  {
1295  const size_t i1( element->index() );
1296  const ET1 v1( element->value() );
1297  ++element;
1298  const size_t i2( element->index() );
1299  const ET1 v2( element->value() );
1300  ++element;
1301  const size_t i3( element->index() );
1302  const ET1 v3( element->value() );
1303  ++element;
1304  const size_t i4( element->index() );
1305  const ET1 v4( element->value() );
1306  ++element;
1307 
1308  BLAZE_INTERNAL_ASSERT( i1 < i2 && i2 < i3 && i3 < i4, "Invalid sparse matrix index detected" );
1309 
1310  const size_t jbegin( ( IsUpper<MT5>::value )
1311  ?( max( jj, ( IsStrictlyUpper<MT5>::value ? i1+1UL : i1 ) ) )
1312  :( jj ) );
1313  const size_t jend( ( IsLower<MT5>::value )
1314  ?( min( jtmp, ( IsStrictlyLower<MT5>::value ? i4 : i4+1UL ) ) )
1315  :( jtmp ) );
1316 
1317  if( IsTriangular<MT5>::value && jbegin >= jend )
1318  continue;
1319 
1320  BLAZE_INTERNAL_ASSERT( jbegin <= jend, "Invalid loop indices detected" );
1321 
1322  const size_t jnum( jend - jbegin );
1323  const size_t jpos( jbegin + ( jnum & size_t(-4) ) );
1324  BLAZE_INTERNAL_ASSERT( ( jbegin + jnum - ( jnum % 4UL ) ) == jpos, "Invalid end calculation" );
1325 
1326  for( size_t j=jbegin; j<jpos; j+=4UL ) {
1327  C(i,j ) -= v1 * B(i1,j ) + v2 * B(i2,j ) + v3 * B(i3,j ) + v4 * B(i4,j );
1328  C(i,j+1UL) -= v1 * B(i1,j+1UL) + v2 * B(i2,j+1UL) + v3 * B(i3,j+1UL) + v4 * B(i4,j+1UL);
1329  C(i,j+2UL) -= v1 * B(i1,j+2UL) + v2 * B(i2,j+2UL) + v3 * B(i3,j+2UL) + v4 * B(i4,j+2UL);
1330  C(i,j+3UL) -= v1 * B(i1,j+3UL) + v2 * B(i2,j+3UL) + v3 * B(i3,j+3UL) + v4 * B(i4,j+3UL);
1331  }
1332  for( size_t j=jpos; j<jend; ++j ) {
1333  C(i,j) -= v1 * B(i1,j) + v2 * B(i2,j) + v3 * B(i3,j) + v4 * B(i4,j);
1334  }
1335  }
1336 
1337  for( ; element!=end; ++element )
1338  {
1339  const size_t i1( element->index() );
1340  const ET1 v1( element->value() );
1341 
1342  const size_t jbegin( ( IsUpper<MT5>::value )
1343  ?( max( jj, ( IsStrictlyUpper<MT5>::value ? i1+1UL : i1 ) ) )
1344  :( jj ) );
1345  const size_t jend( ( IsLower<MT5>::value )
1346  ?( min( jtmp, ( IsStrictlyLower<MT5>::value ? i1 : i1+1UL ) ) )
1347  :( jtmp ) );
1348 
1349  if( IsTriangular<MT5>::value && jbegin >= jend )
1350  continue;
1351 
1352  BLAZE_INTERNAL_ASSERT( jbegin <= jend, "Invalid loop indices detected" );
1353 
1354  const size_t jnum( jend - jbegin );
1355  const size_t jpos( jbegin + ( jnum & size_t(-4) ) );
1356  BLAZE_INTERNAL_ASSERT( ( jbegin + jnum - ( jnum % 4UL ) ) == jpos, "Invalid end calculation" );
1357 
1358  for( size_t j=jbegin; j<jpos; j+=4UL ) {
1359  C(i,j ) -= v1 * B(i1,j );
1360  C(i,j+1UL) -= v1 * B(i1,j+1UL);
1361  C(i,j+2UL) -= v1 * B(i1,j+2UL);
1362  C(i,j+3UL) -= v1 * B(i1,j+3UL);
1363  }
1364  for( size_t j=jpos; j<jend; ++j ) {
1365  C(i,j) -= v1 * B(i1,j);
1366  }
1367  }
1368  }
1369  }
1370  }
1372  //**********************************************************************************************
1373 
1374  //**Vectorized subtraction assignment to dense matrices*****************************************
1388  template< typename MT3 // Type of the left-hand side target matrix
1389  , typename MT4 // Type of the left-hand side matrix operand
1390  , typename MT5 > // Type of the right-hand side matrix operand
1391  static inline EnableIf_< UseVectorizedKernel<MT3,MT4,MT5> >
1392  selectSubAssignKernel( MT3& C, const MT4& A, const MT5& B )
1393  {
1394  typedef ConstIterator_<MT4> ConstIterator;
1395 
1396  const bool remainder( !IsPadded<MT3>::value || !IsPadded<MT5>::value );
1397 
1398  for( size_t i=0UL; i<A.rows(); ++i )
1399  {
1400  const ConstIterator end( A.end(i) );
1401  ConstIterator element( A.begin(i) );
1402 
1403  const size_t nonzeros( A.nonZeros(i) );
1404  const size_t kpos( nonzeros & size_t(-4) );
1405  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
1406 
1407  for( size_t k=0UL; k<kpos; k+=4UL )
1408  {
1409  const size_t i1( element->index() );
1410  const ET1 v1( element->value() );
1411  ++element;
1412  const size_t i2( element->index() );
1413  const ET1 v2( element->value() );
1414  ++element;
1415  const size_t i3( element->index() );
1416  const ET1 v3( element->value() );
1417  ++element;
1418  const size_t i4( element->index() );
1419  const ET1 v4( element->value() );
1420  ++element;
1421 
1422  BLAZE_INTERNAL_ASSERT( i1 < i2 && i2 < i3 && i3 < i4, "Invalid sparse matrix index detected" );
1423 
1424  const SIMDType xmm1( set( v1 ) );
1425  const SIMDType xmm2( set( v2 ) );
1426  const SIMDType xmm3( set( v3 ) );
1427  const SIMDType xmm4( set( v4 ) );
1428 
1429  const size_t jbegin( ( IsUpper<MT5>::value )
1430  ?( ( IsStrictlyUpper<MT5>::value ? i1+1UL : i1 ) & size_t(-SIMDSIZE) )
1431  :( 0UL ) );
1432  const size_t jend( ( IsLower<MT5>::value )
1433  ?( IsStrictlyLower<MT5>::value ? i4 : i4+1UL )
1434  :( B.columns() ) );
1435  BLAZE_INTERNAL_ASSERT( jbegin <= jend, "Invalid loop indices detected" );
1436 
1437  const size_t jpos( remainder ? ( jend & size_t(-SIMDSIZE) ) : jend );
1438  BLAZE_INTERNAL_ASSERT( !remainder || ( jend - ( jend % (SIMDSIZE) ) ) == jpos, "Invalid end calculation" );
1439 
1440  size_t j( jbegin );
1441 
1442  for( ; j<jpos; j+=SIMDSIZE ) {
1443  C.store( i, j, C.load(i,j) - xmm1 * B.load(i1,j) - xmm2 * B.load(i2,j) - xmm3 * B.load(i3,j) - xmm4 * B.load(i4,j) );
1444  }
1445  for( ; remainder && j<jend; ++j ) {
1446  C(i,j) -= v1 * B(i1,j) + v2 * B(i2,j) + v3 * B(i3,j) + v4 * B(i4,j);
1447  }
1448  }
1449 
1450  for( ; element!=end; ++element )
1451  {
1452  const size_t i1( element->index() );
1453  const ET1 v1( element->value() );
1454 
1455  const SIMDType xmm1( set( v1 ) );
1456 
1457  const size_t jbegin( ( IsUpper<MT5>::value )
1458  ?( ( IsStrictlyUpper<MT5>::value ? i1+1UL : i1 ) & size_t(-SIMDSIZE) )
1459  :( 0UL ) );
1460  const size_t jend( ( IsLower<MT5>::value )
1461  ?( IsStrictlyLower<MT5>::value ? i1 : i1+1UL )
1462  :( B.columns() ) );
1463  BLAZE_INTERNAL_ASSERT( jbegin <= jend, "Invalid loop indices detected" );
1464 
1465  const size_t jpos( remainder ? ( jend & size_t(-SIMDSIZE) ) : jend );
1466  BLAZE_INTERNAL_ASSERT( !remainder || ( jend - ( jend % (SIMDSIZE) ) ) == jpos, "Invalid end calculation" );
1467 
1468  size_t j( jbegin );
1469 
1470  for( ; j<jpos; j+=SIMDSIZE ) {
1471  C.store( i, j, C.load(i,j) - xmm1 * B.load(i1,j) );
1472  }
1473  for( ; remainder && j<jend; ++j ) {
1474  C(i,j) -= v1 * B(i1,j);
1475  }
1476  }
1477  }
1478  }
1480  //**********************************************************************************************
1481 
1482  //**Subtraction assignment to sparse matrices***************************************************
1483  // No special implementation for the subtraction assignment to sparse matrices.
1484  //**********************************************************************************************
1485 
1486  //**Multiplication assignment to dense matrices*************************************************
1487  // No special implementation for the multiplication assignment to dense matrices.
1488  //**********************************************************************************************
1489 
1490  //**Multiplication assignment to sparse matrices************************************************
1491  // No special implementation for the multiplication assignment to sparse matrices.
1492  //**********************************************************************************************
1493 
1494  //**SMP assignment to dense matrices************************************************************
1509  template< typename MT // Type of the target dense matrix
1510  , bool SO > // Storage order of the target dense matrix
1511  friend inline EnableIf_< IsEvaluationRequired<MT,MT1,MT2> >
1512  smpAssign( DenseMatrix<MT,SO>& lhs, const SMatDMatMultExpr& rhs )
1513  {
1515 
1516  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1517  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1518 
1519  LT A( rhs.lhs_ ); // Evaluation of the left-hand side sparse matrix operand
1520  RT B( rhs.rhs_ ); // Evaluation of the right-hand side dense matrix operand
1521 
1522  BLAZE_INTERNAL_ASSERT( A.rows() == rhs.lhs_.rows() , "Invalid number of rows" );
1523  BLAZE_INTERNAL_ASSERT( A.columns() == rhs.lhs_.columns(), "Invalid number of columns" );
1524  BLAZE_INTERNAL_ASSERT( B.rows() == rhs.rhs_.rows() , "Invalid number of rows" );
1525  BLAZE_INTERNAL_ASSERT( B.columns() == rhs.rhs_.columns(), "Invalid number of columns" );
1526  BLAZE_INTERNAL_ASSERT( A.rows() == (~lhs).rows() , "Invalid number of rows" );
1527  BLAZE_INTERNAL_ASSERT( B.columns() == (~lhs).columns() , "Invalid number of columns" );
1528 
1529  smpAssign( ~lhs, A * B );
1530  }
1532  //**********************************************************************************************
1533 
1534  //**SMP assignment to sparse matrices***********************************************************
1549  template< typename MT // Type of the target sparse matrix
1550  , bool SO > // Storage order of the target sparse matrix
1551  friend inline EnableIf_< IsEvaluationRequired<MT,MT1,MT2> >
1552  smpAssign( SparseMatrix<MT,SO>& lhs, const SMatDMatMultExpr& rhs )
1553  {
1555 
1556  typedef IfTrue_< SO, OppositeType, ResultType > TmpType;
1557 
1563  BLAZE_CONSTRAINT_MUST_BE_REFERENCE_TYPE( CompositeType_<TmpType> );
1564 
1565  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1566  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1567 
1568  const TmpType tmp( rhs );
1569  smpAssign( ~lhs, tmp );
1570  }
1572  //**********************************************************************************************
1573 
1574  //**SMP addition assignment to dense matrices***************************************************
1590  template< typename MT // Type of the target dense matrix
1591  , bool SO > // Storage order of the target dense matrix
1592  friend inline EnableIf_< IsEvaluationRequired<MT,MT1,MT2> >
1593  smpAddAssign( DenseMatrix<MT,SO>& lhs, const SMatDMatMultExpr& rhs )
1594  {
1596 
1597  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1598  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1599 
1600  LT A( rhs.lhs_ ); // Evaluation of the left-hand side sparse matrix operand
1601  RT B( rhs.rhs_ ); // Evaluation of the right-hand side dense matrix operand
1602 
1603  BLAZE_INTERNAL_ASSERT( A.rows() == rhs.lhs_.rows() , "Invalid number of rows" );
1604  BLAZE_INTERNAL_ASSERT( A.columns() == rhs.lhs_.columns(), "Invalid number of columns" );
1605  BLAZE_INTERNAL_ASSERT( B.rows() == rhs.rhs_.rows() , "Invalid number of rows" );
1606  BLAZE_INTERNAL_ASSERT( B.columns() == rhs.rhs_.columns(), "Invalid number of columns" );
1607  BLAZE_INTERNAL_ASSERT( A.rows() == (~lhs).rows() , "Invalid number of rows" );
1608  BLAZE_INTERNAL_ASSERT( B.columns() == (~lhs).columns() , "Invalid number of columns" );
1609 
1610  smpAddAssign( ~lhs, A * B );
1611  }
1613  //**********************************************************************************************
1614 
1615  //**SMP addition assignment to sparse matrices**************************************************
1616  // No special implementation for the SMP addition assignment to sparse matrices.
1617  //**********************************************************************************************
1618 
1619  //**SMP subtraction assignment to dense matrices************************************************
1635  template< typename MT // Type of the target dense matrix
1636  , bool SO > // Storage order of the target dense matrix
1637  friend inline EnableIf_< IsEvaluationRequired<MT,MT1,MT2> >
1638  smpSubAssign( DenseMatrix<MT,SO>& lhs, const SMatDMatMultExpr& rhs )
1639  {
1641 
1642  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1643  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1644 
1645  LT A( rhs.lhs_ ); // Evaluation of the left-hand side sparse matrix operand
1646  RT B( rhs.rhs_ ); // Evaluation of the right-hand side dense matrix operand
1647 
1648  BLAZE_INTERNAL_ASSERT( A.rows() == rhs.lhs_.rows() , "Invalid number of rows" );
1649  BLAZE_INTERNAL_ASSERT( A.columns() == rhs.lhs_.columns(), "Invalid number of columns" );
1650  BLAZE_INTERNAL_ASSERT( B.rows() == rhs.rhs_.rows() , "Invalid number of rows" );
1651  BLAZE_INTERNAL_ASSERT( B.columns() == rhs.rhs_.columns(), "Invalid number of columns" );
1652  BLAZE_INTERNAL_ASSERT( A.rows() == (~lhs).rows() , "Invalid number of rows" );
1653  BLAZE_INTERNAL_ASSERT( B.columns() == (~lhs).columns() , "Invalid number of columns" );
1654 
1655  smpSubAssign( ~lhs, A * B );
1656  }
1658  //**********************************************************************************************
1659 
1660  //**SMP subtraction assignment to sparse matrices***********************************************
1661  // No special implementation for the SMP subtraction assignment to sparse matrices.
1662  //**********************************************************************************************
1663 
1664  //**SMP multiplication assignment to dense matrices*********************************************
1665  // No special implementation for the SMP multiplication assignment to dense matrices.
1666  //**********************************************************************************************
1667 
1668  //**SMP multiplication assignment to sparse matrices********************************************
1669  // No special implementation for the SMP multiplication assignment to sparse matrices.
1670  //**********************************************************************************************
1671 
1672  //**Compile time checks*************************************************************************
1680  //**********************************************************************************************
1681 };
1682 //*************************************************************************************************
1683 
1684 
1685 
1686 
1687 //=================================================================================================
1688 //
1689 // GLOBAL BINARY ARITHMETIC OPERATORS
1690 //
1691 //=================================================================================================
1692 
1693 //*************************************************************************************************
1722 template< typename T1 // Type of the left-hand side sparse matrix
1723  , typename T2 > // Type of the right-hand side dense matrix
1724 inline const SMatDMatMultExpr<T1,T2>
1726 {
1728 
1729  if( (~lhs).columns() != (~rhs).rows() ) {
1730  BLAZE_THROW_INVALID_ARGUMENT( "Matrix sizes do not match" );
1731  }
1732 
1733  return SMatDMatMultExpr<T1,T2>( ~lhs, ~rhs );
1734 }
1735 //*************************************************************************************************
1736 
1737 
1738 
1739 
1740 //=================================================================================================
1741 //
1742 // ROWS SPECIALIZATIONS
1743 //
1744 //=================================================================================================
1745 
1746 //*************************************************************************************************
1748 template< typename MT1, typename MT2 >
1749 struct Rows< SMatDMatMultExpr<MT1,MT2> > : public Rows<MT1>
1750 {};
1752 //*************************************************************************************************
1753 
1754 
1755 
1756 
1757 //=================================================================================================
1758 //
1759 // COLUMNS SPECIALIZATIONS
1760 //
1761 //=================================================================================================
1762 
1763 //*************************************************************************************************
1765 template< typename MT1, typename MT2 >
1766 struct Columns< SMatDMatMultExpr<MT1,MT2> > : public Columns<MT2>
1767 {};
1769 //*************************************************************************************************
1770 
1771 
1772 
1773 
1774 //=================================================================================================
1775 //
1776 // ISALIGNED SPECIALIZATIONS
1777 //
1778 //=================================================================================================
1779 
1780 //*************************************************************************************************
1782 template< typename MT1, typename MT2 >
1783 struct IsAligned< SMatDMatMultExpr<MT1,MT2> >
1784  : public BoolConstant< IsAligned<MT2>::value >
1785 {};
1787 //*************************************************************************************************
1788 
1789 
1790 
1791 
1792 //=================================================================================================
1793 //
1794 // ISLOWER SPECIALIZATIONS
1795 //
1796 //=================================================================================================
1797 
1798 //*************************************************************************************************
1800 template< typename MT1, typename MT2 >
1801 struct IsLower< SMatDMatMultExpr<MT1,MT2> >
1802  : public BoolConstant< And< IsLower<MT1>, IsLower<MT2> >::value >
1803 {};
1805 //*************************************************************************************************
1806 
1807 
1808 
1809 
1810 //=================================================================================================
1811 //
1812 // ISUNILOWER SPECIALIZATIONS
1813 //
1814 //=================================================================================================
1815 
1816 //*************************************************************************************************
1818 template< typename MT1, typename MT2 >
1819 struct IsUniLower< SMatDMatMultExpr<MT1,MT2> >
1820  : public BoolConstant< And< IsUniLower<MT1>, IsUniLower<MT2> >::value >
1821 {};
1823 //*************************************************************************************************
1824 
1825 
1826 
1827 
1828 //=================================================================================================
1829 //
1830 // ISSTRICTLYLOWER SPECIALIZATIONS
1831 //
1832 //=================================================================================================
1833 
1834 //*************************************************************************************************
1836 template< typename MT1, typename MT2 >
1837 struct IsStrictlyLower< SMatDMatMultExpr<MT1,MT2> >
1838  : public BoolConstant< Or< And< IsStrictlyLower<MT1>, IsLower<MT2> >
1839  , And< IsStrictlyLower<MT2>, IsLower<MT1> > >::value >
1840 {};
1842 //*************************************************************************************************
1843 
1844 
1845 
1846 
1847 //=================================================================================================
1848 //
1849 // ISUPPER SPECIALIZATIONS
1850 //
1851 //=================================================================================================
1852 
1853 //*************************************************************************************************
1855 template< typename MT1, typename MT2 >
1856 struct IsUpper< SMatDMatMultExpr<MT1,MT2> >
1857  : public BoolConstant< And< IsUpper<MT1>, IsUpper<MT2> >::value >
1858 {};
1860 //*************************************************************************************************
1861 
1862 
1863 
1864 
1865 //=================================================================================================
1866 //
1867 // ISUNIUPPER SPECIALIZATIONS
1868 //
1869 //=================================================================================================
1870 
1871 //*************************************************************************************************
1873 template< typename MT1, typename MT2 >
1874 struct IsUniUpper< SMatDMatMultExpr<MT1,MT2> >
1875  : public BoolConstant< And< IsUniUpper<MT1>, IsUniUpper<MT2> >::value >
1876 {};
1878 //*************************************************************************************************
1879 
1880 
1881 
1882 
1883 //=================================================================================================
1884 //
1885 // ISSTRICTLYUPPER SPECIALIZATIONS
1886 //
1887 //=================================================================================================
1888 
1889 //*************************************************************************************************
1891 template< typename MT1, typename MT2 >
1892 struct IsStrictlyUpper< SMatDMatMultExpr<MT1,MT2> >
1893  : public BoolConstant< Or< And< IsStrictlyUpper<MT1>, IsUpper<MT2> >
1894  , And< IsStrictlyUpper<MT2>, IsUpper<MT1> > >::value >
1895 {};
1897 //*************************************************************************************************
1898 
1899 
1900 
1901 
1902 //=================================================================================================
1903 //
1904 // EXPRESSION TRAIT SPECIALIZATIONS
1905 //
1906 //=================================================================================================
1907 
1908 //*************************************************************************************************
1910 template< typename MT1, typename MT2, typename VT >
1911 struct DMatDVecMultExprTrait< SMatDMatMultExpr<MT1,MT2>, VT >
1912 {
1913  public:
1914  //**********************************************************************************************
1915  using Type = If_< And< IsSparseMatrix<MT1>, IsRowMajorMatrix<MT1>
1916  , IsDenseMatrix<MT2>, IsRowMajorMatrix<MT2>
1917  , IsDenseVector<VT>, IsColumnVector<VT> >
1918  , SMatDVecMultExprTrait_< MT1, DMatDVecMultExprTrait_<MT2,VT> >
1919  , INVALID_TYPE >;
1920  //**********************************************************************************************
1921 };
1923 //*************************************************************************************************
1924 
1925 
1926 //*************************************************************************************************
1928 template< typename MT1, typename MT2, typename VT >
1929 struct DMatSVecMultExprTrait< SMatDMatMultExpr<MT1,MT2>, VT >
1930 {
1931  public:
1932  //**********************************************************************************************
1933  using Type = If_< And< IsSparseMatrix<MT1>, IsRowMajorMatrix<MT1>
1934  , IsDenseMatrix<MT2>, IsRowMajorMatrix<MT2>
1935  , IsSparseVector<VT>, IsColumnVector<VT> >
1936  , SMatDVecMultExprTrait_< MT1, DMatSVecMultExprTrait_<MT2,VT> >
1937  , INVALID_TYPE >;
1938  //**********************************************************************************************
1939 };
1941 //*************************************************************************************************
1942 
1943 
1944 //*************************************************************************************************
1946 template< typename VT, typename MT1, typename MT2 >
1947 struct TDVecDMatMultExprTrait< VT, SMatDMatMultExpr<MT1,MT2> >
1948 {
1949  public:
1950  //**********************************************************************************************
1951  using Type = If_< And< IsDenseVector<VT>, IsRowVector<VT>
1952  , IsSparseMatrix<MT1>, IsRowMajorMatrix<MT1>
1953  , IsDenseMatrix<MT2>, IsRowMajorMatrix<MT2> >
1954  , TDVecDMatMultExprTrait_< TDVecSMatMultExprTrait_<VT,MT1>, MT2 >
1955  , INVALID_TYPE >;
1956  //**********************************************************************************************
1957 };
1959 //*************************************************************************************************
1960 
1961 
1962 //*************************************************************************************************
1964 template< typename VT, typename MT1, typename MT2 >
1965 struct TSVecDMatMultExprTrait< VT, SMatDMatMultExpr<MT1,MT2> >
1966 {
1967  public:
1968  //**********************************************************************************************
1969  using Type = If_< And< IsSparseVector<VT>, IsRowVector<VT>
1970  , IsSparseMatrix<MT1>, IsRowMajorMatrix<MT1>
1971  , IsDenseMatrix<MT2>, IsRowMajorMatrix<MT2> >
1972  , TSVecDMatMultExprTrait_< TSVecSMatMultExprTrait_<VT,MT1>, MT2 >
1973  , INVALID_TYPE >;
1974  //**********************************************************************************************
1975 };
1977 //*************************************************************************************************
1978 
1979 
1980 //*************************************************************************************************
1982 template< typename MT1, typename MT2, bool AF >
1983 struct SubmatrixExprTrait< SMatDMatMultExpr<MT1,MT2>, AF >
1984 {
1985  public:
1986  //**********************************************************************************************
1987  using Type = MultExprTrait_< SubmatrixExprTrait_<const MT1,AF>
1988  , SubmatrixExprTrait_<const MT2,AF> >;
1989  //**********************************************************************************************
1990 };
1992 //*************************************************************************************************
1993 
1994 
1995 //*************************************************************************************************
1997 template< typename MT1, typename MT2 >
1998 struct RowExprTrait< SMatDMatMultExpr<MT1,MT2> >
1999 {
2000  public:
2001  //**********************************************************************************************
2002  using Type = MultExprTrait_< RowExprTrait_<const MT1>, MT2 >;
2003  //**********************************************************************************************
2004 };
2006 //*************************************************************************************************
2007 
2008 
2009 //*************************************************************************************************
2011 template< typename MT1, typename MT2 >
2012 struct ColumnExprTrait< SMatDMatMultExpr<MT1,MT2> >
2013 {
2014  public:
2015  //**********************************************************************************************
2016  using Type = MultExprTrait_< MT1, ColumnExprTrait_<const MT2> >;
2017  //**********************************************************************************************
2018 };
2020 //*************************************************************************************************
2021 
2022 } // namespace blaze
2023 
2024 #endif
#define BLAZE_THROW_INVALID_ARGUMENT(MESSAGE)
Macro for the emission of a std::invalid_argument exception.This macro encapsulates the default way o...
Definition: Exception.h:235
Header file for auxiliary alias declarations.
Compile time check whether the given type is a computational expression template.This type trait clas...
Definition: IsComputation.h:72
Header file for mathematical functions.
constexpr bool useOptimizedKernels
Configuration switch for optimized kernels.This configuration switch enables/disables all optimized c...
Definition: Optimizations.h:84
LeftOperand lhs_
Left-hand side sparse matrix of the multiplication expression.
Definition: SMatDMatMultExpr.h:421
Header file for the SMatDVecMultExprTrait class template.
Header file for the Rows type trait.
Header file for the IsUniUpper type trait.
const DMatDMatMultExpr< T1, T2 > operator*(const DenseMatrix< T1, false > &lhs, const DenseMatrix< T2, false > &rhs)
Multiplication operator for the multiplication of two row-major dense matrices ( ).
Definition: DMatDMatMultExpr.h:7800
Compile time check for triangular matrix types.This type trait tests whether or not the given templat...
Definition: IsTriangular.h:87
Header file for basic type definitions.
CompositeType_< MT2 > CT2
Composite type of the right-hand side dense matrix expression.
Definition: SMatDMatMultExpr.h:138
EnableIf_< IsDenseMatrix< MT1 > > smpSubAssign(Matrix< MT1, SO1 > &lhs, const Matrix< MT2, SO2 > &rhs)
Default implementation of the SMP subtraction assignment of a matrix to dense matrix.
Definition: DenseMatrix.h:160
Header file for the IsSparseMatrix type trait.
Header file for the serial shim.
Header file for the IsDiagonal type trait.
BLAZE_ALWAYS_INLINE size_t size(const Vector< VT, TF > &vector) noexcept
Returns the current size/dimension of the vector.
Definition: Vector.h:258
bool isAligned() const noexcept
Returns whether the operands of the expression are properly aligned in memory.
Definition: SMatDMatMultExpr.h:404
#define BLAZE_CONSTRAINT_MUST_BE_DENSE_MATRIX_TYPE(T)
Constraint on the data type.In case the given data type T is not a dense, N-dimensional matrix type...
Definition: DenseMatrix.h:61
Header file for the ColumnExprTrait class template.
ResultType_< MT1 > RT1
Result type of the left-hand side sparse matrix expression.
Definition: SMatDMatMultExpr.h:133
IfTrue_< evaluateLeft, const RT1, CT1 > LT
Type for the assignment of the left-hand side sparse matrix operand.
Definition: SMatDMatMultExpr.h:232
If_< IsExpression< MT1 >, const MT1, const MT1 & > LeftOperand
Composite type of the left-hand side sparse matrix expression.
Definition: SMatDMatMultExpr.h:226
BLAZE_ALWAYS_INLINE MT::Iterator begin(Matrix< MT, SO > &matrix, size_t i)
Returns an iterator to the first element of row/column i.
Definition: Matrix.h:188
Availability of a SIMD multiplication for the given data types.Depending on the available instruction...
Definition: HasSIMDMult.h:162
typename SIMDTrait< T >::Type SIMDTrait_
Auxiliary alias declaration for the SIMDTrait class template.The SIMDTrait_ alias declaration provide...
Definition: SIMDTrait.h:315
void reset(const DiagonalProxy< MT > &proxy)
Resetting the represented element to the default initial values.
Definition: DiagonalProxy.h:533
Header file for the IsRowVector type trait.
Header file for the And class template.
const ElementType_< MT > min(const DenseMatrix< MT, SO > &dm)
Returns the smallest element of the dense matrix.
Definition: DenseMatrix.h:1669
Compile time check for lower triangular matrices.This type trait tests whether or not the given templ...
Definition: IsLower.h:88
Availability of a SIMD addition for the given data types.Depending on the available instruction set (...
Definition: HasSIMDAdd.h:162
Header file for the TDVecSMatMultExprTrait class template.
const DMatSerialExpr< MT, SO > serial(const DenseMatrix< MT, SO > &dm)
Forces the serial evaluation of the given dense matrix expression dm.
Definition: DMatSerialExpr.h:723
typename MultTrait< T1, T2 >::Type MultTrait_
Auxiliary alias declaration for the MultTrait class template.The MultTrait_ alias declaration provide...
Definition: MultTrait.h:245
Header file for the Computation base class.
Header file for the MatMatMultExpr base class.
Expression object for sparse matrix-dense matrix multiplications.The SMatDMatMultExpr class represent...
Definition: Forward.h:90
Compile time check for upper triangular matrices.This type trait tests whether or not the given templ...
Definition: IsUpper.h:88
If_< IsExpression< MT2 >, const MT2, const MT2 & > RightOperand
Composite type of the right-hand side dense matrix expression.
Definition: SMatDMatMultExpr.h:229
Constraints on the storage order of matrix types.
Header file for the RequiresEvaluation type trait.
ElementType_< ResultType > ElementType
Resulting element type.
Definition: SMatDMatMultExpr.h:220
System settings for performance optimizations.
Header file for the TSVecSMatMultExprTrait class template.
Header file for the IsUniLower type trait.
typename T::ResultType ResultType_
Alias declaration for nested ResultType type definitions.The ResultType_ alias declaration provides a...
Definition: Aliases.h:323
const ElementType_< MT > max(const DenseMatrix< MT, SO > &dm)
Returns the largest element of the dense matrix.
Definition: DenseMatrix.h:1716
EnableIf_< IsDenseMatrix< MT1 > > smpAddAssign(Matrix< MT1, SO1 > &lhs, const Matrix< MT2, SO2 > &rhs)
Default implementation of the SMP addition assignment of a matrix to a dense matrix.
Definition: DenseMatrix.h:129
DisableIf_< Or< IsComputation< MT >, IsTransExpr< MT > >, ColumnExprTrait_< MT > > column(Matrix< MT, SO > &matrix, size_t index)
Creating a view on a specific column of the given matrix.
Definition: Column.h:126
Base class for dense matrices.The DenseMatrix class is a base class for all dense matrix classes...
Definition: DenseMatrix.h:70
Base class for sparse matrices.The SparseMatrix class is a base class for all sparse matrix classes...
Definition: Forward.h:109
bool isDefault(const DiagonalProxy< MT > &proxy)
Returns whether the represented element is in default state.
Definition: DiagonalProxy.h:573
typename IfTrue< Condition, T1, T2 >::Type IfTrue_
Auxiliary alias declaration for the IfTrue class template.The IfTrue_ alias declaration provides a co...
Definition: If.h:109
Constraint on the data type.
Constraint on the data type.
Header file for the MultExprTrait class template.
SIMDTrait_< ElementType > SIMDType
Resulting SIMD element type.
Definition: SMatDMatMultExpr.h:221
Compile time check to query the requirement to evaluate an expression.Via this type trait it is possi...
Definition: RequiresEvaluation.h:72
SubvectorExprTrait_< VT, unaligned > subvector(Vector< VT, TF > &vector, size_t index, size_t size)
Creating a view on a specific subvector of the given vector.
Definition: Subvector.h:152
bool isAliased(const T *alias) const noexcept
Returns whether the expression is aliased with the given address alias.
Definition: SMatDMatMultExpr.h:394
typename T::CompositeType CompositeType_
Alias declaration for nested CompositeType type definitions.The CompositeType_ alias declaration prov...
Definition: Aliases.h:83
Header file for the multiplication trait.
Header file for the IsStrictlyUpper type trait.
Namespace of the Blaze C++ math library.
Definition: Blaze.h:57
Header file for the If class template.
SMatDMatMultExpr(const MT1 &lhs, const MT2 &rhs) noexcept
Constructor for the SMatDMatMultExpr class.
Definition: SMatDMatMultExpr.h:261
Compile time check for row-major matrix types.This type trait tests whether or not the given template...
Definition: IsRowMajorMatrix.h:83
#define BLAZE_CONSTRAINT_MUST_BE_COLUMN_MAJOR_MATRIX_TYPE(T)
Constraint on the data type.In case the given data type T is not a column-major dense or sparse matri...
Definition: ColumnMajorMatrix.h:61
ReturnType operator()(size_t i, size_t j) const
2D-access to the matrix elements.
Definition: SMatDMatMultExpr.h:276
bool canAlias(const T *alias) const noexcept
Returns whether the expression can alias with the given address alias.
Definition: SMatDMatMultExpr.h:382
const Element * ConstIterator
Iterator over constant elements.
Definition: CompressedMatrix.h:2647
EnableIf_< IsDenseMatrix< MT1 > > smpAssign(Matrix< MT1, SO1 > &lhs, const Matrix< MT2, SO2 > &rhs)
Default implementation of the SMP assignment of a matrix to a dense matrix.
Definition: DenseMatrix.h:98
Header file for the Or class template.
ElementType_< RT2 > ET2
Element type of the right-hand side dense matrix expression.
Definition: SMatDMatMultExpr.h:136
#define BLAZE_THROW_OUT_OF_RANGE(MESSAGE)
Macro for the emission of a std::out_of_range exception.This macro encapsulates the default way of Bl...
Definition: Exception.h:331
Header file for the HasSIMDAdd type trait.
Header file for the DenseMatrix base class.
Header file for the Columns type trait.
typename T::ElementType ElementType_
Alias declaration for nested ElementType type definitions.The ElementType_ alias declaration provides...
Definition: Aliases.h:163
Header file for all SIMD functionality.
MultTrait_< RT1, RT2 > ResultType
Result type for expression template evaluations.
Definition: SMatDMatMultExpr.h:217
OppositeType_< ResultType > OppositeType
Result type with opposite storage order for expression template evaluations.
Definition: SMatDMatMultExpr.h:218
Header file for the DMatDVecMultExprTrait class template.
LeftOperand leftOperand() const noexcept
Returns the left-hand side sparse matrix operand.
Definition: SMatDMatMultExpr.h:360
Header file for the IsLower type trait.
Header file for the IsAligned type trait.
Compile time check for diagonal matrices.This type trait tests whether or not the given template para...
Definition: IsDiagonal.h:90
#define BLAZE_CONSTRAINT_MUST_BE_REFERENCE_TYPE(T)
Constraint on the data type.In case the given data type T is not a reference type, a compilation error is created.
Definition: Reference.h:60
Header file for the IsTriangular type trait.
ResultType_< MT2 > RT2
Result type of the right-hand side dense matrix expression.
Definition: SMatDMatMultExpr.h:134
Constraints on the storage order of matrix types.
Header file for the exception macros of the math module.
Compile time check for strictly upper triangular matrices.This type trait tests whether or not the gi...
Definition: IsStrictlyUpper.h:86
BLAZE_ALWAYS_INLINE MT::Iterator end(Matrix< MT, SO > &matrix, size_t i)
Returns an iterator just past the last element of row/column i.
Definition: Matrix.h:254
Header file for the RowExprTrait class template.
Header file for all forward declarations for expression class templates.
Header file for the IsDenseMatrix type trait.
TransposeType_< ResultType > TransposeType
Transpose type for expression template evaluations.
Definition: SMatDMatMultExpr.h:219
Header file for the EnableIf class template.
Header file for the IsStrictlyLower type trait.
#define BLAZE_CONSTRAINT_MUST_FORM_VALID_MATMATMULTEXPR(T1, T2)
Constraint on the data type.In case the given data types T1 and T2 do not form a valid matrix/matrix ...
Definition: MatMatMultExpr.h:109
BLAZE_ALWAYS_INLINE const EnableIf_< And< IsIntegral< T >, HasSize< T, 1UL > >, If_< IsSigned< T >, SIMDint8, SIMDuint8 > > set(T value) noexcept
Sets all values in the vector to the given 1-byte integral value.
Definition: Set.h:76
DisableIf_< Or< IsComputation< MT >, IsTransExpr< MT > >, RowExprTrait_< MT > > row(Matrix< MT, SO > &matrix, size_t index)
Creating a view on a specific row of the given matrix.
Definition: Row.h:126
ReturnType at(size_t i, size_t j) const
Checked access to the matrix elements.
Definition: SMatDMatMultExpr.h:324
Header file for the IsSparseVector type trait.
Header file for the SubmatrixExprTrait class template.
#define BLAZE_CONSTRAINT_MUST_BE_ROW_MAJOR_MATRIX_TYPE(T)
Constraint on the data type.In case the given data type T is not a row-major dense or sparse matrix t...
Definition: RowMajorMatrix.h:61
Header file for the HasSIMDMult type trait.
const ElementType ReturnType
Return type for expression template evaluations.
Definition: SMatDMatMultExpr.h:222
Header file for run time assertion macros.
Utility type for generic codes.
RightOperand rightOperand() const noexcept
Returns the right-hand side dense matrix operand.
Definition: SMatDMatMultExpr.h:370
typename If< T1, T2, T3 >::Type If_
Auxiliary alias declaration for the If class template.The If_ alias declaration provides a convenient...
Definition: If.h:160
Header file for the reset shim.
SMatDMatMultExpr< MT1, MT2 > This
Type of this SMatDMatMultExpr instance.
Definition: SMatDMatMultExpr.h:216
Constraints on the storage order of matrix types.
RightOperand rhs_
Right-hand side dense matrix of the multiplication expression.
Definition: SMatDMatMultExpr.h:422
IfTrue_< evaluateRight, const RT2, CT2 > RT
Type for the assignment of the right-hand side dense matrix operand.
Definition: SMatDMatMultExpr.h:235
IntegralConstant< bool, B > BoolConstant
Generic wrapper for a compile time constant boolean value.The BoolConstant class template represents ...
Definition: IntegralConstant.h:100
Header file for the RemoveReference type trait.
typename EnableIf< Condition, T >::Type EnableIf_
Auxiliary alias declaration for the EnableIf class template.The EnableIf_ alias declaration provides ...
Definition: EnableIf.h:223
typename T::OppositeType OppositeType_
Alias declaration for nested OppositeType type definitions.The OppositeType_ alias declaration provid...
Definition: Aliases.h:243
#define BLAZE_CONSTRAINT_MATRICES_MUST_HAVE_SAME_STORAGE_ORDER(T1, T2)
Constraint on the data type.In case either of the two given data types T1 or T2 is not a matrix type ...
Definition: StorageOrder.h:84
Header file for the IsDenseVector type trait.
ElementType_< RT1 > ET1
Element type of the left-hand side sparse matrix expression.
Definition: SMatDMatMultExpr.h:135
Compile time check for strictly lower triangular matrices.This type trait tests whether or not the gi...
Definition: IsStrictlyLower.h:86
typename T::ConstIterator ConstIterator_
Alias declaration for nested ConstIterator type definitions.The ConstIterator_ alias declaration prov...
Definition: Aliases.h:103
Header file for the AreSIMDCombinable type trait.
Header file for the IsRowMajorMatrix type trait.
size_t rows() const noexcept
Returns the current number of rows of the matrix.
Definition: SMatDMatMultExpr.h:340
Header file for the IsComputation type trait class.
Header file for the TDVecDMatMultExprTrait class template.
Compile time logical or evaluation.The Or class template performs at compile time a logical or ('&&')...
Definition: Or.h:101
#define BLAZE_FUNCTION_TRACE
Function trace macro.This macro can be used to reliably trace function calls. In case function tracin...
Definition: FunctionTrace.h:157
Header file for the IntegralConstant class template.
Header file for the TSVecDMatMultExprTrait class template.
CompositeType_< MT1 > CT1
Composite type of the left-hand side sparse matrix expression.
Definition: SMatDMatMultExpr.h:137
const ResultType CompositeType
Data type for composite expression templates.
Definition: SMatDMatMultExpr.h:223
typename T::TransposeType TransposeType_
Alias declaration for nested TransposeType type definitions.The TransposeType_ alias declaration prov...
Definition: Aliases.h:403
Header file for the IsUpper type trait.
Header file for the DMatSVecMultExprTrait class template.
bool canSMPAssign() const noexcept
Returns whether the expression can be used in SMP assignments.
Definition: SMatDMatMultExpr.h:414
Header file for the IsColumnVector type trait.
Constraint on the data type.
Header file for the IsResizable type trait.
Header file for the thresholds for matrix/vector and matrix/matrix multiplications.
#define BLAZE_INTERNAL_ASSERT(expr, msg)
Run time assertion macro for internal checks.In case of an invalid run time expression, the program execution is terminated. The BLAZE_INTERNAL_ASSERT macro can be disabled by setting the BLAZE_USER_ASSERTION flag to zero or by defining NDEBUG during the compilation.
Definition: Assert.h:101
#define BLAZE_CONSTRAINT_MUST_BE_SPARSE_MATRIX_TYPE(T)
Constraint on the data type.In case the given data type T is not a sparse, N-dimensional matrix type...
Definition: SparseMatrix.h:61
size_t columns() const noexcept
Returns the current number of columns of the matrix.
Definition: SMatDMatMultExpr.h:350
Header file for the IsExpression type trait class.
Header file for the FunctionTrace class.