SMatTDMatMultExpr.h
Go to the documentation of this file.
1 //=================================================================================================
33 //=================================================================================================
34 
35 #ifndef _BLAZE_MATH_EXPRESSIONS_SMATTDMATMULTEXPR_H_
36 #define _BLAZE_MATH_EXPRESSIONS_SMATTDMATMULTEXPR_H_
37 
38 
39 //*************************************************************************************************
40 // Includes
41 //*************************************************************************************************
42 
43 #include <blaze/math/Aliases.h>
52 #include <blaze/math/Exception.h>
64 #include <blaze/math/shims/Reset.h>
87 #include <blaze/math/views/Check.h>
92 #include <blaze/util/Assert.h>
93 #include <blaze/util/DisableIf.h>
94 #include <blaze/util/EnableIf.h>
97 #include <blaze/util/mpl/If.h>
98 #include <blaze/util/TrueType.h>
99 #include <blaze/util/Types.h>
101 
102 
103 namespace blaze {
104 
105 //=================================================================================================
106 //
107 // CLASS SMATTDMATMULTEXPR
108 //
109 //=================================================================================================
110 
111 //*************************************************************************************************
118 template< typename MT1 // Type of the left-hand side sparse matrix
119  , typename MT2 // Type of the right-hand side dense matrix
120  , bool SF // Symmetry flag
121  , bool HF // Hermitian flag
122  , bool LF // Lower flag
123  , bool UF > // Upper flag
124 class SMatTDMatMultExpr
125  : public MatMatMultExpr< DenseMatrix< SMatTDMatMultExpr<MT1,MT2,SF,HF,LF,UF>, true > >
126  , private Computation
127 {
128  private:
129  //**Type definitions****************************************************************************
136  //**********************************************************************************************
137 
138  //**********************************************************************************************
140  static constexpr bool evaluateLeft = ( IsComputation_v<MT1> || RequiresEvaluation_v<MT1> );
141  //**********************************************************************************************
142 
143  //**********************************************************************************************
145  static constexpr bool evaluateRight = ( IsComputation_v<MT2> || RequiresEvaluation_v<MT2> );
146  //**********************************************************************************************
147 
148  //**********************************************************************************************
149  static constexpr bool SYM = ( SF && !( HF || LF || UF ) );
150  static constexpr bool HERM = ( HF && !( LF || UF ) );
151  static constexpr bool LOW = ( LF || ( ( SF || HF ) && UF ) );
152  static constexpr bool UPP = ( UF || ( ( SF || HF ) && LF ) );
153  //**********************************************************************************************
154 
155  //**********************************************************************************************
157 
162  template< typename T1, typename T2, typename T3 >
163  static constexpr bool CanExploitSymmetry_v = IsSymmetric_v<T3>;
165  //**********************************************************************************************
166 
167  //**********************************************************************************************
169 
173  template< typename T1, typename T2, typename T3 >
174  static constexpr bool IsEvaluationRequired_v =
175  ( ( evaluateLeft || evaluateRight ) && CanExploitSymmetry_v<T1,T2,T3> );
177  //**********************************************************************************************
178 
179  //**********************************************************************************************
181 
184  template< typename T1, typename T2, typename T3 >
185  static constexpr bool UseOptimizedKernel_v =
186  ( useOptimizedKernels &&
187  !IsDiagonal_v<T3> &&
188  !IsResizable_v< ElementType_t<T1> > &&
189  !IsResizable_v<ET1> );
191  //**********************************************************************************************
192 
193  //**********************************************************************************************
195 
198  using ForwardFunctor = If_t< HERM
199  , DeclHerm
200  , If_t< SYM
201  , DeclSym
202  , If_t< LOW
203  , If_t< UPP
204  , DeclDiag
205  , DeclLow >
206  , If_t< UPP
207  , DeclUpp
208  , Noop > > > >;
210  //**********************************************************************************************
211 
212  public:
213  //**Type definitions****************************************************************************
216 
219 
221  using ResultType = typename If_t< HERM
223  , If_t< SYM
225  , If_t< LOW
226  , If_t< UPP
229  , If_t< UPP
231  , MultTrait<RT1,RT2> > > > >::Type;
232 
236  using ReturnType = const ElementType;
237  using CompositeType = const ResultType;
238 
240  using LeftOperand = If_t< IsExpression_v<MT1>, const MT1, const MT1& >;
241 
243  using RightOperand = If_t< IsExpression_v<MT2>, const MT2, const MT2& >;
244 
247 
250  //**********************************************************************************************
251 
252  //**Compilation flags***************************************************************************
254  static constexpr bool simdEnabled = false;
255 
257  static constexpr bool smpAssignable =
259  //**********************************************************************************************
260 
261  //**Constructor*********************************************************************************
267  explicit inline SMatTDMatMultExpr( const MT1& lhs, const MT2& rhs ) noexcept
268  : lhs_( lhs ) // Left-hand side sparse matrix of the multiplication expression
269  , rhs_( rhs ) // Right-hand side dense matrix of the multiplication expression
270  {
271  BLAZE_INTERNAL_ASSERT( lhs.columns() == rhs.rows(), "Invalid matrix sizes" );
272  }
273  //**********************************************************************************************
274 
275  //**Access operator*****************************************************************************
282  inline ReturnType operator()( size_t i, size_t j ) const {
283  BLAZE_INTERNAL_ASSERT( i < lhs_.rows() , "Invalid row access index" );
284  BLAZE_INTERNAL_ASSERT( j < rhs_.columns(), "Invalid column access index" );
285 
286  if( IsDiagonal_v<MT1> ) {
287  return lhs_(i,i) * rhs_(i,j);
288  }
289  else if( IsDiagonal_v<MT2> ) {
290  return lhs_(i,j) * rhs_(j,j);
291  }
292  else if( IsTriangular_v<MT1> || IsTriangular_v<MT2> ) {
293  const size_t begin( ( IsUpper_v<MT1> )
294  ?( ( IsLower_v<MT2> )
295  ?( max( ( IsStrictlyUpper_v<MT1> ? i+1UL : i )
296  , ( IsStrictlyLower_v<MT2> ? j+1UL : j ) ) )
297  :( IsStrictlyUpper_v<MT1> ? i+1UL : i ) )
298  :( ( IsLower_v<MT2> )
299  ?( IsStrictlyLower_v<MT2> ? j+1UL : j )
300  :( 0UL ) ) );
301  const size_t end( ( IsLower_v<MT1> )
302  ?( ( IsUpper_v<MT2> )
303  ?( min( ( IsStrictlyLower_v<MT1> ? i : i+1UL )
304  , ( IsStrictlyUpper_v<MT2> ? j : j+1UL ) ) )
305  :( IsStrictlyLower_v<MT1> ? i : i+1UL ) )
306  :( ( IsUpper_v<MT2> )
307  ?( IsStrictlyUpper_v<MT2> ? j : j+1UL )
308  :( lhs_.columns() ) ) );
309 
310  if( begin >= end ) return ElementType();
311 
312  const size_t n( end - begin );
313 
314  return subvector( row( lhs_, i, unchecked ), begin, n, unchecked ) *
315  subvector( column( rhs_, j, unchecked ), begin, n, unchecked );
316  }
317  else {
318  return row( lhs_, i, unchecked ) * column( rhs_, j, unchecked );
319  }
320  }
321  //**********************************************************************************************
322 
323  //**At function*********************************************************************************
331  inline ReturnType at( size_t i, size_t j ) const {
332  if( i >= lhs_.rows() ) {
333  BLAZE_THROW_OUT_OF_RANGE( "Invalid row access index" );
334  }
335  if( j >= rhs_.columns() ) {
336  BLAZE_THROW_OUT_OF_RANGE( "Invalid column access index" );
337  }
338  return (*this)(i,j);
339  }
340  //**********************************************************************************************
341 
342  //**Rows function*******************************************************************************
347  inline size_t rows() const noexcept {
348  return lhs_.rows();
349  }
350  //**********************************************************************************************
351 
352  //**Columns function****************************************************************************
357  inline size_t columns() const noexcept {
358  return rhs_.columns();
359  }
360  //**********************************************************************************************
361 
362  //**Left operand access*************************************************************************
367  inline LeftOperand leftOperand() const noexcept {
368  return lhs_;
369  }
370  //**********************************************************************************************
371 
372  //**Right operand access************************************************************************
377  inline RightOperand rightOperand() const noexcept {
378  return rhs_;
379  }
380  //**********************************************************************************************
381 
382  //**********************************************************************************************
388  template< typename T >
389  inline bool canAlias( const T* alias ) const noexcept {
390  return ( lhs_.isAliased( alias ) || rhs_.isAliased( alias ) );
391  }
392  //**********************************************************************************************
393 
394  //**********************************************************************************************
400  template< typename T >
401  inline bool isAliased( const T* alias ) const noexcept {
402  return ( lhs_.isAliased( alias ) || rhs_.isAliased( alias ) );
403  }
404  //**********************************************************************************************
405 
406  //**********************************************************************************************
411  inline bool isAligned() const noexcept {
412  return rhs_.isAligned();
413  }
414  //**********************************************************************************************
415 
416  //**********************************************************************************************
421  inline bool canSMPAssign() const noexcept {
422  return ( rows() * columns() >= SMP_SMATTDMATMULT_THRESHOLD ) && !IsDiagonal_v<MT2>;
423  }
424  //**********************************************************************************************
425 
426  private:
427  //**Member variables****************************************************************************
430  //**********************************************************************************************
431 
432  //**Assignment to dense matrices****************************************************************
445  template< typename MT // Type of the target dense matrix
446  , bool SO > // Storage order of the target dense matrix
447  friend inline auto assign( DenseMatrix<MT,SO>& lhs, const SMatTDMatMultExpr& rhs )
449  {
451 
452  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
453  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
454 
455  LT A( serial( rhs.lhs_ ) ); // Evaluation of the left-hand side sparse matrix operand
456  RT B( serial( rhs.rhs_ ) ); // Evaluation of the right-hand side dense matrix operand
457 
458  BLAZE_INTERNAL_ASSERT( A.rows() == rhs.lhs_.rows() , "Invalid number of rows" );
459  BLAZE_INTERNAL_ASSERT( A.columns() == rhs.lhs_.columns(), "Invalid number of columns" );
460  BLAZE_INTERNAL_ASSERT( B.rows() == rhs.rhs_.rows() , "Invalid number of rows" );
461  BLAZE_INTERNAL_ASSERT( B.columns() == rhs.rhs_.columns(), "Invalid number of columns" );
462  BLAZE_INTERNAL_ASSERT( A.rows() == (~lhs).rows() , "Invalid number of rows" );
463  BLAZE_INTERNAL_ASSERT( B.columns() == (~lhs).columns() , "Invalid number of columns" );
464 
465  SMatTDMatMultExpr::selectAssignKernel( ~lhs, A, B );
466  }
468  //**********************************************************************************************
469 
470  //**Default assignment to dense matrices********************************************************
484  template< typename MT3 // Type of the left-hand side target matrix
485  , typename MT4 // Type of the left-hand side matrix operand
486  , typename MT5 > // Type of the right-hand side matrix operand
487  static inline auto selectAssignKernel( MT3& C, const MT4& A, const MT5& B )
489  {
490  const size_t M( A.rows() );
491  const size_t N( B.columns() );
492 
493  BLAZE_INTERNAL_ASSERT( !( SYM || HERM || LOW || UPP ) || M == N, "Broken invariant detected" );
494 
495  {
496  size_t j( 0UL );
497 
498  for( ; (j+4UL) <= N; j+=4UL ) {
499  for( size_t i=( SYM || HERM || LOW ? j : 0UL ); i<( UPP ? j+4UL : M ); ++i )
500  {
501  const auto end( ( IsUpper_v<MT5> )
502  ?( IsStrictlyUpper_v<MT5> ? A.lowerBound(i,j+4UL) : A.upperBound(i,j+4UL) )
503  :( A.end(i) ) );
504  auto element( ( IsLower_v<MT5> )
505  ?( IsStrictlyLower_v<MT5> ? A.upperBound(i,j) : A.lowerBound(i,j) )
506  :( A.begin(i) ) );
507 
508  if( element == end ) {
509  reset( C(i,j ) );
510  reset( C(i,j+1UL) );
511  reset( C(i,j+2UL) );
512  reset( C(i,j+3UL) );
513  continue;
514  }
515 
516  C(i,j ) = element->value() * B(element->index(),j );
517  C(i,j+1UL) = element->value() * B(element->index(),j+1UL);
518  C(i,j+2UL) = element->value() * B(element->index(),j+2UL);
519  C(i,j+3UL) = element->value() * B(element->index(),j+3UL);
520  ++element;
521  for( ; element!=end; ++element ) {
522  C(i,j ) += element->value() * B(element->index(),j );
523  C(i,j+1UL) += element->value() * B(element->index(),j+1UL);
524  C(i,j+2UL) += element->value() * B(element->index(),j+2UL);
525  C(i,j+3UL) += element->value() * B(element->index(),j+3UL);
526  }
527  }
528  }
529 
530  for( ; (j+2UL) <= N; j+=2UL ) {
531  for( size_t i=( SYM || HERM || LOW ? j : 0UL ); i<( UPP ? j+2UL : M ); ++i )
532  {
533  const auto end( ( IsUpper_v<MT5> )
534  ?( IsStrictlyUpper_v<MT5> ? A.lowerBound(i,j+2UL) : A.upperBound(i,j+2UL) )
535  :( A.end(i) ) );
536  auto element( ( IsLower_v<MT5> )
537  ?( IsStrictlyLower_v<MT5> ? A.upperBound(i,j) : A.lowerBound(i,j) )
538  :( A.begin(i) ) );
539 
540  if( element == end ) {
541  reset( C(i,j ) );
542  reset( C(i,j+1UL) );
543  continue;
544  }
545 
546  C(i,j ) = element->value() * B(element->index(),j );
547  C(i,j+1UL) = element->value() * B(element->index(),j+1UL);
548  ++element;
549  for( ; element!=end; ++element ) {
550  C(i,j ) += element->value() * B(element->index(),j );
551  C(i,j+1UL) += element->value() * B(element->index(),j+1UL);
552  }
553  }
554  }
555 
556  for( ; j<N; ++j ) {
557  for( size_t i=( SYM || HERM || LOW ? j : 0UL ); i<( UPP ? j+1UL : M ); ++i )
558  {
559  const auto end( ( IsUpper_v<MT5> )
560  ?( IsStrictlyUpper_v<MT5> ? A.lowerBound(i,j) : A.upperBound(i,j) )
561  :( A.end(i) ) );
562  auto element( ( IsLower_v<MT5> )
563  ?( IsStrictlyLower_v<MT5> ? A.upperBound(i,j) : A.lowerBound(i,j) )
564  :( A.begin(i) ) );
565 
566  if( element == end ) {
567  reset( C(i,j) );
568  continue;
569  }
570 
571  C(i,j) = element->value() * B(element->index(),j);
572  ++element;
573  for( ; element!=end; ++element ) {
574  C(i,j) += element->value() * B(element->index(),j);
575  }
576  }
577  }
578  }
579 
580  if( SYM || HERM ) {
581  for( size_t j=1UL; j<N; ++j ) {
582  for( size_t i=0UL; i<j; ++i ) {
583  C(i,j) = HERM ? conj( C(j,i) ) : C(j,i);
584  }
585  }
586  }
587  else if( LOW && !UPP ) {
588  for( size_t j=1UL; j<N; ++j ) {
589  for( size_t i=0UL; i<j; ++i ) {
590  reset( C(i,j) );
591  }
592  }
593  }
594  else if( !LOW && UPP ) {
595  for( size_t i=1UL; i<M; ++i ) {
596  for( size_t j=0UL; j<i; ++j ) {
597  reset( C(i,j) );
598  }
599  }
600  }
601  }
603  //**********************************************************************************************
604 
605  //**Optimized assignment to dense matrices******************************************************
619  template< typename MT3 // Type of the left-hand side target matrix
620  , typename MT4 // Type of the left-hand side matrix operand
621  , typename MT5 > // Type of the right-hand side matrix operand
622  static inline auto selectAssignKernel( MT3& C, const MT4& A, const MT5& B )
623  -> EnableIf_t< UseOptimizedKernel_v<MT3,MT4,MT5> >
624  {
625  const size_t M( A.rows() );
626  const size_t N( B.columns() );
627 
628  BLAZE_INTERNAL_ASSERT( !( SYM || HERM || LOW || UPP ) || M == N, "Broken invariant detected" );
629 
630  reset( C );
631 
632  {
633  size_t j( 0UL );
634 
635  for( ; (j+4UL) <= N; j+=4UL ) {
636  for( size_t i=( SYM || HERM || LOW ? j : 0UL ); i<( UPP ? j+4UL : M ); ++i )
637  {
638  const auto end( ( IsUpper_v<MT5> )
639  ?( IsStrictlyUpper_v<MT5> ? A.lowerBound(i,j+4UL) : A.upperBound(i,j+4UL) )
640  :( A.end(i) ) );
641  auto element( ( IsLower_v<MT5> )
642  ?( IsStrictlyLower_v<MT5> ? A.upperBound(i,j) : A.lowerBound(i,j) )
643  :( A.begin(i) ) );
644 
645  const size_t nonzeros( end - element );
646  const size_t kpos( nonzeros & size_t(-4) );
647  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
648 
649  for( size_t k=0UL; k<kpos; k+=4UL )
650  {
651  const size_t i1( element->index() );
652  const ET1 v1( element->value() );
653  ++element;
654  const size_t i2( element->index() );
655  const ET1 v2( element->value() );
656  ++element;
657  const size_t i3( element->index() );
658  const ET1 v3( element->value() );
659  ++element;
660  const size_t i4( element->index() );
661  const ET1 v4( element->value() );
662  ++element;
663 
664  BLAZE_INTERNAL_ASSERT( i1 < i2 && i2 < i3 && i3 < i4, "Invalid sparse matrix index detected" );
665 
666  C(i,j ) += v1 * B(i1,j ) + v2 * B(i2,j ) + v3 * B(i3,j ) + v4 * B(i4,j );
667  C(i,j+1UL) += v1 * B(i1,j+1UL) + v2 * B(i2,j+1UL) + v3 * B(i3,j+1UL) + v4 * B(i4,j+1UL);
668  C(i,j+2UL) += v1 * B(i1,j+2UL) + v2 * B(i2,j+2UL) + v3 * B(i3,j+2UL) + v4 * B(i4,j+2UL);
669  C(i,j+3UL) += v1 * B(i1,j+3UL) + v2 * B(i2,j+3UL) + v3 * B(i3,j+3UL) + v4 * B(i4,j+3UL);
670  }
671 
672  for( ; element!=end; ++element )
673  {
674  const size_t i1( element->index() );
675  const ET1 v1( element->value() );
676 
677  C(i,j ) += v1 * B(i1,j );
678  C(i,j+1UL) += v1 * B(i1,j+1UL);
679  C(i,j+2UL) += v1 * B(i1,j+2UL);
680  C(i,j+3UL) += v1 * B(i1,j+3UL);
681  }
682  }
683  }
684 
685  for( ; (j+2UL) <= N; j+=2UL ) {
686  for( size_t i=( SYM || HERM || LOW ? j : 0UL ); i<( UPP ? j+2UL : M ); ++i )
687  {
688  const auto end( ( IsUpper_v<MT5> )
689  ?( IsStrictlyUpper_v<MT5> ? A.lowerBound(i,j+2UL) : A.upperBound(i,j+2UL) )
690  :( A.end(i) ) );
691  auto element( ( IsLower_v<MT5> )
692  ?( IsStrictlyLower_v<MT5> ? A.upperBound(i,j) : A.lowerBound(i,j) )
693  :( A.begin(i) ) );
694 
695  const size_t nonzeros( end - element );
696  const size_t kpos( nonzeros & size_t(-4) );
697  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
698 
699  for( size_t k=0UL; k<kpos; k+=4UL )
700  {
701  const size_t i1( element->index() );
702  const ET1 v1( element->value() );
703  ++element;
704  const size_t i2( element->index() );
705  const ET1 v2( element->value() );
706  ++element;
707  const size_t i3( element->index() );
708  const ET1 v3( element->value() );
709  ++element;
710  const size_t i4( element->index() );
711  const ET1 v4( element->value() );
712  ++element;
713 
714  BLAZE_INTERNAL_ASSERT( i1 < i2 && i2 < i3 && i3 < i4, "Invalid sparse matrix index detected" );
715 
716  C(i,j ) += v1 * B(i1,j ) + v2 * B(i2,j ) + v3 * B(i3,j ) + v4 * B(i4,j );
717  C(i,j+1UL) += v1 * B(i1,j+1UL) + v2 * B(i2,j+1UL) + v3 * B(i3,j+1UL) + v4 * B(i4,j+1UL);
718  }
719 
720  for( ; element!=end; ++element )
721  {
722  const size_t i1( element->index() );
723  const ET1 v1( element->value() );
724 
725  C(i,j ) += v1 * B(i1,j );
726  C(i,j+1UL) += v1 * B(i1,j+1UL);
727  }
728  }
729  }
730 
731  for( ; j<N; ++j ) {
732  for( size_t i=( SYM || HERM || LOW ? j : 0UL ); i<( UPP ? j+1UL : M ); ++i )
733  {
734  const auto end( ( IsUpper_v<MT5> )
735  ?( IsStrictlyUpper_v<MT5> ? A.lowerBound(i,j) : A.upperBound(i,j) )
736  :( A.end(i) ) );
737  auto element( ( IsLower_v<MT5> )
738  ?( IsStrictlyLower_v<MT5> ? A.upperBound(i,j) : A.lowerBound(i,j) )
739  :( A.begin(i) ) );
740 
741  const size_t nonzeros( end - element );
742  const size_t kpos( nonzeros & size_t(-4) );
743  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
744 
745  for( size_t k=0UL; k<kpos; k+=4UL )
746  {
747  const size_t i1( element->index() );
748  const ET1 v1( element->value() );
749  ++element;
750  const size_t i2( element->index() );
751  const ET1 v2( element->value() );
752  ++element;
753  const size_t i3( element->index() );
754  const ET1 v3( element->value() );
755  ++element;
756  const size_t i4( element->index() );
757  const ET1 v4( element->value() );
758  ++element;
759 
760  BLAZE_INTERNAL_ASSERT( i1 < i2 && i2 < i3 && i3 < i4, "Invalid sparse matrix index detected" );
761 
762  C(i,j) += v1 * B(i1,j) + v2 * B(i2,j) + v3 * B(i3,j) + v4 * B(i4,j);
763  }
764 
765  for( ; element!=end; ++element )
766  {
767  const size_t i1( element->index() );
768  const ET1 v1( element->value() );
769 
770  C(i,j) += v1 * B(i1,j);
771  }
772  }
773  }
774  }
775 
776  if( SYM || HERM ) {
777  for( size_t j=1UL; j<N; ++j ) {
778  for( size_t i=0UL; i<j; ++i ) {
779  C(i,j) = HERM ? conj( C(j,i) ) : C(j,i);
780  }
781  }
782  }
783  }
785  //**********************************************************************************************
786 
787  //**Assignment to sparse matrices***************************************************************
800  template< typename MT // Type of the target sparse matrix
801  , bool SO > // Storage order of the target sparse matrix
802  friend inline auto assign( SparseMatrix<MT,SO>& lhs, const SMatTDMatMultExpr& rhs )
803  -> DisableIf_t< CanExploitSymmetry_v<MT,MT1,MT2> >
804  {
806 
807  using TmpType = If_t< SO, ResultType, OppositeType >;
808 
815 
816  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
817  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
818 
819  const ForwardFunctor fwd;
820 
821  const TmpType tmp( serial( rhs ) );
822  assign( ~lhs, fwd( tmp ) );
823  }
825  //**********************************************************************************************
826 
827  //**Restructuring assignment********************************************************************
842  template< typename MT // Type of the target matrix
843  , bool SO > // Storage order of the target matrix
844  friend inline auto assign( Matrix<MT,SO>& lhs, const SMatTDMatMultExpr& rhs )
845  -> EnableIf_t< CanExploitSymmetry_v<MT,MT1,MT2> >
846  {
848 
849  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
850  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
851 
852  const ForwardFunctor fwd;
853 
854  assign( ~lhs, fwd( rhs.lhs_ * trans( rhs.rhs_ ) ) );
855  }
857  //**********************************************************************************************
858 
859  //**Addition assignment to dense matrices*******************************************************
872  template< typename MT // Type of the target dense matrix
873  , bool SO > // Storage order of the target dense matrix
874  friend inline auto addAssign( DenseMatrix<MT,SO>& lhs, const SMatTDMatMultExpr& rhs )
875  -> DisableIf_t< CanExploitSymmetry_v<MT,MT1,MT2> >
876  {
878 
879  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
880  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
881 
882  LT A( serial( rhs.lhs_ ) ); // Evaluation of the left-hand side sparse matrix operand
883  RT B( serial( rhs.rhs_ ) ); // Evaluation of the right-hand side dense matrix operand
884 
885  BLAZE_INTERNAL_ASSERT( A.rows() == rhs.lhs_.rows() , "Invalid number of rows" );
886  BLAZE_INTERNAL_ASSERT( A.columns() == rhs.lhs_.columns(), "Invalid number of columns" );
887  BLAZE_INTERNAL_ASSERT( B.rows() == rhs.rhs_.rows() , "Invalid number of rows" );
888  BLAZE_INTERNAL_ASSERT( B.columns() == rhs.rhs_.columns(), "Invalid number of columns" );
889  BLAZE_INTERNAL_ASSERT( A.rows() == (~lhs).rows() , "Invalid number of rows" );
890  BLAZE_INTERNAL_ASSERT( B.columns() == (~lhs).columns() , "Invalid number of columns" );
891 
892  SMatTDMatMultExpr::selectAddAssignKernel( ~lhs, A, B );
893  }
895  //**********************************************************************************************
896 
897  //**Default addition assignment to dense matrices***********************************************
911  template< typename MT3 // Type of the left-hand side target matrix
912  , typename MT4 // Type of the left-hand side matrix operand
913  , typename MT5 > // Type of the right-hand side matrix operand
914  static inline auto selectAddAssignKernel( MT3& C, const MT4& A, const MT5& B )
915  -> DisableIf_t< UseOptimizedKernel_v<MT3,MT4,MT5> >
916  {
917  const size_t M( A.rows() );
918  const size_t N( B.columns() );
919 
920  BLAZE_INTERNAL_ASSERT( !( LOW || UPP ) || M == N, "Broken invariant detected" );
921 
922  {
923  size_t j( 0UL );
924 
925  for( ; (j+4UL) <= N; j+=4UL ) {
926  for( size_t i=( LOW ? j : 0UL ); i<( UPP ? j+4UL : M ); ++i )
927  {
928  const auto end( ( IsUpper_v<MT5> )
929  ?( IsStrictlyUpper_v<MT5> ? A.lowerBound(i,j+4UL) : A.upperBound(i,j+4UL) )
930  :( A.end(i) ) );
931  auto element( ( IsLower_v<MT5> )
932  ?( IsStrictlyLower_v<MT5> ? A.upperBound(i,j) : A.lowerBound(i,j) )
933  :( A.begin(i) ) );
934 
935  for( ; element!=end; ++element ) {
936  C(i,j ) += element->value() * B(element->index(),j );
937  C(i,j+1UL) += element->value() * B(element->index(),j+1UL);
938  C(i,j+2UL) += element->value() * B(element->index(),j+2UL);
939  C(i,j+3UL) += element->value() * B(element->index(),j+3UL);
940  }
941  }
942  }
943 
944  for( ; (j+2UL) <= N; j+=2UL ) {
945  for( size_t i=( LOW ? j : 0UL ); i<( UPP ? j+2UL : M ); ++i )
946  {
947  const auto end( ( IsUpper_v<MT5> )
948  ?( IsStrictlyUpper_v<MT5> ? A.lowerBound(i,j+2UL) : A.upperBound(i,j+2UL) )
949  :( A.end(i) ) );
950  auto element( ( IsLower_v<MT5> )
951  ?( IsStrictlyLower_v<MT5> ? A.upperBound(i,j) : A.lowerBound(i,j) )
952  :( A.begin(i) ) );
953 
954  for( ; element!=end; ++element ) {
955  C(i,j ) += element->value() * B(element->index(),j );
956  C(i,j+1UL) += element->value() * B(element->index(),j+1UL);
957  }
958  }
959  }
960 
961  for( ; j<N; ++j ) {
962  for( size_t i=( LOW ? j : 0UL ); i<( UPP ? j+1UL : M ); ++i )
963  {
964  const auto end( ( IsUpper_v<MT5> )
965  ?( IsStrictlyUpper_v<MT5> ? A.lowerBound(i,j) : A.upperBound(i,j) )
966  :( A.end(i) ) );
967  auto element( ( IsLower_v<MT5> )
968  ?( IsStrictlyLower_v<MT5> ? A.upperBound(i,j) : A.lowerBound(i,j) )
969  :( A.begin(i) ) );
970 
971  for( ; element!=end; ++element ) {
972  C(i,j) += element->value() * B(element->index(),j);
973  }
974  }
975  }
976  }
977  }
979  //**********************************************************************************************
980 
981  //**Optimized addition assignment to dense matrices*********************************************
995  template< typename MT3 // Type of the left-hand side target matrix
996  , typename MT4 // Type of the left-hand side matrix operand
997  , typename MT5 > // Type of the right-hand side matrix operand
998  static inline auto selectAddAssignKernel( MT3& C, const MT4& A, const MT5& B )
999  -> EnableIf_t< UseOptimizedKernel_v<MT3,MT4,MT5> >
1000  {
1001  const size_t M( A.rows() );
1002  const size_t N( B.columns() );
1003 
1004  BLAZE_INTERNAL_ASSERT( !( LOW || UPP ) || M == N, "Broken invariant detected" );
1005 
1006  {
1007  size_t j( 0UL );
1008 
1009  for( ; (j+4UL) <= N; j+=4UL ) {
1010  for( size_t i=( LOW ? j : 0UL ); i<( UPP ? j+4UL : M ); ++i )
1011  {
1012  const auto end( ( IsUpper_v<MT5> )
1013  ?( IsStrictlyUpper_v<MT5> ? A.lowerBound(i,j+4UL) : A.upperBound(i,j+4UL) )
1014  :( A.end(i) ) );
1015  auto element( ( IsLower_v<MT5> )
1016  ?( IsStrictlyLower_v<MT5> ? A.upperBound(i,j) : A.lowerBound(i,j) )
1017  :( A.begin(i) ) );
1018 
1019  const size_t nonzeros( end - element );
1020  const size_t kpos( nonzeros & size_t(-4) );
1021  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
1022 
1023  for( size_t k=0UL; k<kpos; k+=4UL )
1024  {
1025  const size_t i1( element->index() );
1026  const ET1 v1( element->value() );
1027  ++element;
1028  const size_t i2( element->index() );
1029  const ET1 v2( element->value() );
1030  ++element;
1031  const size_t i3( element->index() );
1032  const ET1 v3( element->value() );
1033  ++element;
1034  const size_t i4( element->index() );
1035  const ET1 v4( element->value() );
1036  ++element;
1037 
1038  BLAZE_INTERNAL_ASSERT( i1 < i2 && i2 < i3 && i3 < i4, "Invalid sparse matrix index detected" );
1039 
1040  C(i,j ) += v1 * B(i1,j ) + v2 * B(i2,j ) + v3 * B(i3,j ) + v4 * B(i4,j );
1041  C(i,j+1UL) += v1 * B(i1,j+1UL) + v2 * B(i2,j+1UL) + v3 * B(i3,j+1UL) + v4 * B(i4,j+1UL);
1042  C(i,j+2UL) += v1 * B(i1,j+2UL) + v2 * B(i2,j+2UL) + v3 * B(i3,j+2UL) + v4 * B(i4,j+2UL);
1043  C(i,j+3UL) += v1 * B(i1,j+3UL) + v2 * B(i2,j+3UL) + v3 * B(i3,j+3UL) + v4 * B(i4,j+3UL);
1044  }
1045 
1046  for( ; element!=end; ++element )
1047  {
1048  const size_t i1( element->index() );
1049  const ET1 v1( element->value() );
1050 
1051  C(i,j ) += v1 * B(i1,j );
1052  C(i,j+1UL) += v1 * B(i1,j+1UL);
1053  C(i,j+2UL) += v1 * B(i1,j+2UL);
1054  C(i,j+3UL) += v1 * B(i1,j+3UL);
1055  }
1056  }
1057  }
1058 
1059  for( ; (j+2UL) <= N; j+=2UL ) {
1060  for( size_t i=( LOW ? j : 0UL ); i<( UPP ? j+2UL : M ); ++i )
1061  {
1062  const auto end( ( IsUpper_v<MT5> )
1063  ?( IsStrictlyUpper_v<MT5> ? A.lowerBound(i,j+2UL) : A.upperBound(i,j+2UL) )
1064  :( A.end(i) ) );
1065  auto element( ( IsLower_v<MT5> )
1066  ?( IsStrictlyLower_v<MT5> ? A.upperBound(i,j) : A.lowerBound(i,j) )
1067  :( A.begin(i) ) );
1068 
1069  const size_t nonzeros( end - element );
1070  const size_t kpos( nonzeros & size_t(-4) );
1071  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
1072 
1073  for( size_t k=0UL; k<kpos; k+=4UL )
1074  {
1075  const size_t i1( element->index() );
1076  const ET1 v1( element->value() );
1077  ++element;
1078  const size_t i2( element->index() );
1079  const ET1 v2( element->value() );
1080  ++element;
1081  const size_t i3( element->index() );
1082  const ET1 v3( element->value() );
1083  ++element;
1084  const size_t i4( element->index() );
1085  const ET1 v4( element->value() );
1086  ++element;
1087 
1088  BLAZE_INTERNAL_ASSERT( i1 < i2 && i2 < i3 && i3 < i4, "Invalid sparse matrix index detected" );
1089 
1090  C(i,j ) += v1 * B(i1,j ) + v2 * B(i2,j ) + v3 * B(i3,j ) + v4 * B(i4,j );
1091  C(i,j+1UL) += v1 * B(i1,j+1UL) + v2 * B(i2,j+1UL) + v3 * B(i3,j+1UL) + v4 * B(i4,j+1UL);
1092  }
1093 
1094  for( ; element!=end; ++element )
1095  {
1096  const size_t i1( element->index() );
1097  const ET1 v1( element->value() );
1098 
1099  C(i,j ) += v1 * B(i1,j );
1100  C(i,j+1UL) += v1 * B(i1,j+1UL);
1101  }
1102  }
1103  }
1104 
1105  for( ; j<N; ++j ) {
1106  for( size_t i=( LOW ? j : 0UL ); i<( UPP ? j+1UL : M ); ++i )
1107  {
1108  const auto end( ( IsUpper_v<MT5> )
1109  ?( IsStrictlyUpper_v<MT5> ? A.lowerBound(i,j) : A.upperBound(i,j) )
1110  :( A.end(i) ) );
1111  auto element( ( IsLower_v<MT5> )
1112  ?( IsStrictlyLower_v<MT5> ? A.upperBound(i,j) : A.lowerBound(i,j) )
1113  :( A.begin(i) ) );
1114 
1115  const size_t nonzeros( end - element );
1116  const size_t kpos( nonzeros & size_t(-4) );
1117  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
1118 
1119  for( size_t k=0UL; k<kpos; k+=4UL )
1120  {
1121  const size_t i1( element->index() );
1122  const ET1 v1( element->value() );
1123  ++element;
1124  const size_t i2( element->index() );
1125  const ET1 v2( element->value() );
1126  ++element;
1127  const size_t i3( element->index() );
1128  const ET1 v3( element->value() );
1129  ++element;
1130  const size_t i4( element->index() );
1131  const ET1 v4( element->value() );
1132  ++element;
1133 
1134  BLAZE_INTERNAL_ASSERT( i1 < i2 && i2 < i3 && i3 < i4, "Invalid sparse matrix index detected" );
1135 
1136  C(i,j) += v1 * B(i1,j) + v2 * B(i2,j) + v3 * B(i3,j) + v4 * B(i4,j);
1137  }
1138 
1139  for( ; element!=end; ++element )
1140  {
1141  const size_t i1( element->index() );
1142  const ET1 v1( element->value() );
1143 
1144  C(i,j) += v1 * B(i1,j);
1145  }
1146  }
1147  }
1148  }
1149  }
1151  //**********************************************************************************************
1152 
1153  //**Restructuring addition assignment***********************************************************
1168  template< typename MT // Type of the target matrix
1169  , bool SO > // Storage order of the target matrix
1170  friend inline auto addAssign( Matrix<MT,SO>& lhs, const SMatTDMatMultExpr& rhs )
1171  -> EnableIf_t< CanExploitSymmetry_v<MT,MT1,MT2> >
1172  {
1174 
1175  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1176  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1177 
1178  const ForwardFunctor fwd;
1179 
1180  addAssign( ~lhs, fwd( rhs.lhs_ * trans( rhs.rhs_ ) ) );
1181  }
1183  //**********************************************************************************************
1184 
1185  //**Addition assignment to sparse matrices******************************************************
1186  // No special implementation for the addition assignment to sparse matrices.
1187  //**********************************************************************************************
1188 
1189  //**Subtraction assignment to dense matrices****************************************************
1202  template< typename MT // Type of the target dense matrix
1203  , bool SO > // Storage order of the target dense matrix
1204  friend inline auto subAssign( DenseMatrix<MT,SO>& lhs, const SMatTDMatMultExpr& rhs )
1205  -> DisableIf_t< CanExploitSymmetry_v<MT,MT1,MT2> >
1206  {
1208 
1209  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1210  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1211 
1212  LT A( serial( rhs.lhs_ ) ); // Evaluation of the left-hand side sparse matrix operand
1213  RT B( serial( rhs.rhs_ ) ); // Evaluation of the right-hand side dense matrix operand
1214 
1215  BLAZE_INTERNAL_ASSERT( A.rows() == rhs.lhs_.rows() , "Invalid number of rows" );
1216  BLAZE_INTERNAL_ASSERT( A.columns() == rhs.lhs_.columns(), "Invalid number of columns" );
1217  BLAZE_INTERNAL_ASSERT( B.rows() == rhs.rhs_.rows() , "Invalid number of rows" );
1218  BLAZE_INTERNAL_ASSERT( B.columns() == rhs.rhs_.columns(), "Invalid number of columns" );
1219  BLAZE_INTERNAL_ASSERT( A.rows() == (~lhs).rows() , "Invalid number of rows" );
1220  BLAZE_INTERNAL_ASSERT( B.columns() == (~lhs).columns() , "Invalid number of columns" );
1221 
1222  SMatTDMatMultExpr::selectSubAssignKernel( ~lhs, A, B );
1223  }
1225  //**********************************************************************************************
1226 
1227  //**Default subtraction assignment to dense matrices********************************************
1241  template< typename MT3 // Type of the left-hand side target matrix
1242  , typename MT4 // Type of the left-hand side matrix operand
1243  , typename MT5 > // Type of the right-hand side matrix operand
1244  static inline auto selectSubAssignKernel( MT3& C, const MT4& A, const MT5& B )
1245  -> DisableIf_t< UseOptimizedKernel_v<MT3,MT4,MT5> >
1246  {
1247  const size_t M( A.rows() );
1248  const size_t N( B.columns() );
1249 
1250  BLAZE_INTERNAL_ASSERT( !( LOW || UPP ) || M == N, "Broken invariant detected" );
1251 
1252  {
1253  size_t j( 0UL );
1254 
1255  for( ; (j+4UL) <= N; j+=4UL ) {
1256  for( size_t i=( LOW ? j : 0UL ); i<( UPP ? j+4UL : M ); ++i )
1257  {
1258  const auto end( ( IsUpper_v<MT5> )
1259  ?( IsStrictlyUpper_v<MT5> ? A.lowerBound(i,j+4UL) : A.upperBound(i,j+4UL) )
1260  :( A.end(i) ) );
1261  auto element( ( IsLower_v<MT5> )
1262  ?( IsStrictlyLower_v<MT5> ? A.upperBound(i,j) : A.lowerBound(i,j) )
1263  :( A.begin(i) ) );
1264 
1265  for( ; element!=end; ++element ) {
1266  C(i,j ) -= element->value() * B(element->index(),j );
1267  C(i,j+1UL) -= element->value() * B(element->index(),j+1UL);
1268  C(i,j+2UL) -= element->value() * B(element->index(),j+2UL);
1269  C(i,j+3UL) -= element->value() * B(element->index(),j+3UL);
1270  }
1271  }
1272  }
1273 
1274  for( ; (j+2UL) <= N; j+=2UL ) {
1275  for( size_t i=( LOW ? j : 0UL ); i<( UPP ? j+2UL : M ); ++i )
1276  {
1277  const auto end( ( IsUpper_v<MT5> )
1278  ?( IsStrictlyUpper_v<MT5> ? A.lowerBound(i,j+2UL) : A.upperBound(i,j+2UL) )
1279  :( A.end(i) ) );
1280  auto element( ( IsLower_v<MT5> )
1281  ?( IsStrictlyLower_v<MT5> ? A.upperBound(i,j) : A.lowerBound(i,j) )
1282  :( A.begin(i) ) );
1283 
1284  for( ; element!=end; ++element ) {
1285  C(i,j ) -= element->value() * B(element->index(),j );
1286  C(i,j+1UL) -= element->value() * B(element->index(),j+1UL);
1287  }
1288  }
1289  }
1290 
1291  for( ; j<N; ++j ) {
1292  for( size_t i=( LOW ? j : 0UL ); i<( UPP ? j+1UL : M ); ++i )
1293  {
1294  const auto end( ( IsUpper_v<MT5> )
1295  ?( IsStrictlyUpper_v<MT5> ? A.lowerBound(i,j) : A.upperBound(i,j) )
1296  :( A.end(i) ) );
1297  auto element( ( IsLower_v<MT5> )
1298  ?( IsStrictlyLower_v<MT5> ? A.upperBound(i,j) : A.lowerBound(i,j) )
1299  :( A.begin(i) ) );
1300 
1301  for( ; element!=end; ++element ) {
1302  C(i,j) -= element->value() * B(element->index(),j);
1303  }
1304  }
1305  }
1306  }
1307  }
1309  //**********************************************************************************************
1310 
1311  //**Optimized subtraction assignment to dense matrices******************************************
1325  template< typename MT3 // Type of the left-hand side target matrix
1326  , typename MT4 // Type of the left-hand side matrix operand
1327  , typename MT5 > // Type of the right-hand side matrix operand
1328  static inline auto selectSubAssignKernel( MT3& C, const MT4& A, const MT5& B )
1329  -> EnableIf_t< UseOptimizedKernel_v<MT3,MT4,MT5> >
1330  {
1331  const size_t M( A.rows() );
1332  const size_t N( B.columns() );
1333 
1334  BLAZE_INTERNAL_ASSERT( !( LOW || UPP ) || M == N, "Broken invariant detected" );
1335 
1336  {
1337  size_t j( 0UL );
1338 
1339  for( ; (j+4UL) <= N; j+=4UL ) {
1340  for( size_t i=( LOW ? j : 0UL ); i<( UPP ? j+4UL : M ); ++i )
1341  {
1342  const auto end( ( IsUpper_v<MT5> )
1343  ?( IsStrictlyUpper_v<MT5> ? A.lowerBound(i,j+4UL) : A.upperBound(i,j+4UL) )
1344  :( A.end(i) ) );
1345  auto element( ( IsLower_v<MT5> )
1346  ?( IsStrictlyLower_v<MT5> ? A.upperBound(i,j) : A.lowerBound(i,j) )
1347  :( A.begin(i) ) );
1348 
1349  const size_t nonzeros( end - element );
1350  const size_t kpos( nonzeros & size_t(-4) );
1351  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
1352 
1353  for( size_t k=0UL; k<kpos; k+=4UL )
1354  {
1355  const size_t i1( element->index() );
1356  const ET1 v1( element->value() );
1357  ++element;
1358  const size_t i2( element->index() );
1359  const ET1 v2( element->value() );
1360  ++element;
1361  const size_t i3( element->index() );
1362  const ET1 v3( element->value() );
1363  ++element;
1364  const size_t i4( element->index() );
1365  const ET1 v4( element->value() );
1366  ++element;
1367 
1368  BLAZE_INTERNAL_ASSERT( i1 < i2 && i2 < i3 && i3 < i4, "Invalid sparse matrix index detected" );
1369 
1370  C(i,j ) -= v1 * B(i1,j ) + v2 * B(i2,j ) + v3 * B(i3,j ) + v4 * B(i4,j );
1371  C(i,j+1UL) -= v1 * B(i1,j+1UL) + v2 * B(i2,j+1UL) + v3 * B(i3,j+1UL) + v4 * B(i4,j+1UL);
1372  C(i,j+2UL) -= v1 * B(i1,j+2UL) + v2 * B(i2,j+2UL) + v3 * B(i3,j+2UL) + v4 * B(i4,j+2UL);
1373  C(i,j+3UL) -= v1 * B(i1,j+3UL) + v2 * B(i2,j+3UL) + v3 * B(i3,j+3UL) + v4 * B(i4,j+3UL);
1374  }
1375 
1376  for( ; element!=end; ++element )
1377  {
1378  const size_t i1( element->index() );
1379  const ET1 v1( element->value() );
1380 
1381  C(i,j ) -= v1 * B(i1,j );
1382  C(i,j+1UL) -= v1 * B(i1,j+1UL);
1383  C(i,j+2UL) -= v1 * B(i1,j+2UL);
1384  C(i,j+3UL) -= v1 * B(i1,j+3UL);
1385  }
1386  }
1387  }
1388 
1389  for( ; (j+2UL) <= N; j+=2UL ) {
1390  for( size_t i=( LOW ? j : 0UL ); i<( UPP ? j+2UL : M ); ++i )
1391  {
1392  const auto end( ( IsUpper_v<MT5> )
1393  ?( IsStrictlyUpper_v<MT5> ? A.lowerBound(i,j+2UL) : A.upperBound(i,j+2UL) )
1394  :( A.end(i) ) );
1395  auto element( ( IsLower_v<MT5> )
1396  ?( IsStrictlyLower_v<MT5> ? A.upperBound(i,j) : A.lowerBound(i,j) )
1397  :( A.begin(i) ) );
1398 
1399  const size_t nonzeros( end - element );
1400  const size_t kpos( nonzeros & size_t(-4) );
1401  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
1402 
1403  for( size_t k=0UL; k<kpos; k+=4UL )
1404  {
1405  const size_t i1( element->index() );
1406  const ET1 v1( element->value() );
1407  ++element;
1408  const size_t i2( element->index() );
1409  const ET1 v2( element->value() );
1410  ++element;
1411  const size_t i3( element->index() );
1412  const ET1 v3( element->value() );
1413  ++element;
1414  const size_t i4( element->index() );
1415  const ET1 v4( element->value() );
1416  ++element;
1417 
1418  BLAZE_INTERNAL_ASSERT( i1 < i2 && i2 < i3 && i3 < i4, "Invalid sparse matrix index detected" );
1419 
1420  C(i,j ) -= v1 * B(i1,j ) + v2 * B(i2,j ) + v3 * B(i3,j ) + v4 * B(i4,j );
1421  C(i,j+1UL) -= v1 * B(i1,j+1UL) + v2 * B(i2,j+1UL) + v3 * B(i3,j+1UL) + v4 * B(i4,j+1UL);
1422  }
1423 
1424  for( ; element!=end; ++element )
1425  {
1426  const size_t i1( element->index() );
1427  const ET1 v1( element->value() );
1428 
1429  C(i,j ) -= v1 * B(i1,j );
1430  C(i,j+1UL) -= v1 * B(i1,j+1UL);
1431  }
1432  }
1433  }
1434 
1435  for( ; j<N; ++j ) {
1436  for( size_t i=( LOW ? j : 0UL ); i<( UPP ? j+1UL : M ); ++i )
1437  {
1438  const auto end( ( IsUpper_v<MT5> )
1439  ?( IsStrictlyUpper_v<MT5> ? A.lowerBound(i,j) : A.upperBound(i,j) )
1440  :( A.end(i) ) );
1441  auto element( ( IsLower_v<MT5> )
1442  ?( IsStrictlyLower_v<MT5> ? A.upperBound(i,j) : A.lowerBound(i,j) )
1443  :( A.begin(i) ) );
1444 
1445  const size_t nonzeros( end - element );
1446  const size_t kpos( nonzeros & size_t(-4) );
1447  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
1448 
1449  for( size_t k=0UL; k<kpos; k+=4UL )
1450  {
1451  const size_t i1( element->index() );
1452  const ET1 v1( element->value() );
1453  ++element;
1454  const size_t i2( element->index() );
1455  const ET1 v2( element->value() );
1456  ++element;
1457  const size_t i3( element->index() );
1458  const ET1 v3( element->value() );
1459  ++element;
1460  const size_t i4( element->index() );
1461  const ET1 v4( element->value() );
1462  ++element;
1463 
1464  BLAZE_INTERNAL_ASSERT( i1 < i2 && i2 < i3 && i3 < i4, "Invalid sparse matrix index detected" );
1465 
1466  C(i,j) -= v1 * B(i1,j) + v2 * B(i2,j) + v3 * B(i3,j) + v4 * B(i4,j);
1467  }
1468 
1469  for( ; element!=end; ++element )
1470  {
1471  const size_t i1( element->index() );
1472  const ET1 v1( element->value() );
1473 
1474  C(i,j) -= v1 * B(i1,j);
1475  }
1476  }
1477  }
1478  }
1479  }
1481  //**********************************************************************************************
1482 
1483  //**Restructuring subtraction assignment********************************************************
1498  template< typename MT // Type of the target matrix
1499  , bool SO > // Storage order of the target matrix
1500  friend inline auto subAssign( Matrix<MT,SO>& lhs, const SMatTDMatMultExpr& rhs )
1501  -> EnableIf_t< CanExploitSymmetry_v<MT,MT1,MT2> >
1502  {
1504 
1505  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1506  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1507 
1508  const ForwardFunctor fwd;
1509 
1510  subAssign( ~lhs, fwd( rhs.lhs_ * trans( rhs.rhs_ ) ) );
1511  }
1513  //**********************************************************************************************
1514 
1515  //**Subtraction assignment to sparse matrices***************************************************
1516  // No special implementation for the subtraction assignment to sparse matrices.
1517  //**********************************************************************************************
1518 
1519  //**Schur product assignment to dense matrices**************************************************
1532  template< typename MT // Type of the target dense matrix
1533  , bool SO > // Storage order of the target dense matrix
1534  friend inline void schurAssign( DenseMatrix<MT,SO>& lhs, const SMatTDMatMultExpr& rhs )
1535  {
1537 
1541 
1542  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1543  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1544 
1545  const ResultType tmp( serial( rhs ) );
1546  schurAssign( ~lhs, tmp );
1547  }
1549  //**********************************************************************************************
1550 
1551  //**Schur product assignment to sparse matrices*************************************************
1552  // No special implementation for the Schur product assignment to sparse matrices.
1553  //**********************************************************************************************
1554 
1555  //**Multiplication assignment to dense matrices*************************************************
1556  // No special implementation for the multiplication assignment to dense matrices.
1557  //**********************************************************************************************
1558 
1559  //**Multiplication assignment to sparse matrices************************************************
1560  // No special implementation for the multiplication assignment to sparse matrices.
1561  //**********************************************************************************************
1562 
1563  //**SMP assignment to dense matrices************************************************************
1578  template< typename MT // Type of the target dense matrix
1579  , bool SO > // Storage order of the target dense matrix
1580  friend inline auto smpAssign( DenseMatrix<MT,SO>& lhs, const SMatTDMatMultExpr& rhs )
1581  -> EnableIf_t< IsEvaluationRequired_v<MT,MT1,MT2> >
1582  {
1584 
1585  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1586  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1587 
1588  LT A( rhs.lhs_ ); // Evaluation of the left-hand side sparse matrix operand
1589  RT B( rhs.rhs_ ); // Evaluation of the right-hand side dense matrix operand
1590 
1591  BLAZE_INTERNAL_ASSERT( A.rows() == rhs.lhs_.rows() , "Invalid number of rows" );
1592  BLAZE_INTERNAL_ASSERT( A.columns() == rhs.lhs_.columns(), "Invalid number of columns" );
1593  BLAZE_INTERNAL_ASSERT( B.rows() == rhs.rhs_.rows() , "Invalid number of rows" );
1594  BLAZE_INTERNAL_ASSERT( B.columns() == rhs.rhs_.columns(), "Invalid number of columns" );
1595  BLAZE_INTERNAL_ASSERT( A.rows() == (~lhs).rows() , "Invalid number of rows" );
1596  BLAZE_INTERNAL_ASSERT( B.columns() == (~lhs).columns() , "Invalid number of columns" );
1597 
1598  smpAssign( ~lhs, A * B );
1599  }
1601  //**********************************************************************************************
1602 
1603  //**SMP assignment to sparse matrices***********************************************************
1618  template< typename MT // Type of the target sparse matrix
1619  , bool SO > // Storage order of the target sparse matrix
1620  friend inline auto smpAssign( SparseMatrix<MT,SO>& lhs, const SMatTDMatMultExpr& rhs )
1621  -> EnableIf_t< IsEvaluationRequired_v<MT,MT1,MT2> >
1622  {
1624 
1625  using TmpType = If_t< SO, ResultType, OppositeType >;
1626 
1633 
1634  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1635  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1636 
1637  const ForwardFunctor fwd;
1638 
1639  const TmpType tmp( rhs );
1640  smpAssign( ~lhs, fwd( tmp ) );
1641  }
1643  //**********************************************************************************************
1644 
1645  //**Restructuring SMP assignment****************************************************************
1660  template< typename MT // Type of the target matrix
1661  , bool SO > // Storage order of the target matrix
1662  friend inline auto smpAssign( Matrix<MT,SO>& lhs, const SMatTDMatMultExpr& rhs )
1663  -> EnableIf_t< CanExploitSymmetry_v<MT,MT1,MT2> >
1664  {
1666 
1667  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1668  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1669 
1670  const ForwardFunctor fwd;
1671 
1672  smpAssign( ~lhs, fwd( rhs.lhs_ * trans( rhs.rhs_ ) ) );
1673  }
1675  //**********************************************************************************************
1676 
1677  //**SMP addition assignment to dense matrices***************************************************
1693  template< typename MT // Type of the target dense matrix
1694  , bool SO > // Storage order of the target dense matrix
1695  friend inline auto smpAddAssign( DenseMatrix<MT,SO>& lhs, const SMatTDMatMultExpr& rhs )
1696  -> EnableIf_t< IsEvaluationRequired_v<MT,MT1,MT2> >
1697  {
1699 
1700  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1701  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1702 
1703  LT A( rhs.lhs_ ); // Evaluation of the left-hand side sparse matrix operand
1704  RT B( rhs.rhs_ ); // Evaluation of the right-hand side dense matrix operand
1705 
1706  BLAZE_INTERNAL_ASSERT( A.rows() == rhs.lhs_.rows() , "Invalid number of rows" );
1707  BLAZE_INTERNAL_ASSERT( A.columns() == rhs.lhs_.columns(), "Invalid number of columns" );
1708  BLAZE_INTERNAL_ASSERT( B.rows() == rhs.rhs_.rows() , "Invalid number of rows" );
1709  BLAZE_INTERNAL_ASSERT( B.columns() == rhs.rhs_.columns(), "Invalid number of columns" );
1710  BLAZE_INTERNAL_ASSERT( A.rows() == (~lhs).rows() , "Invalid number of rows" );
1711  BLAZE_INTERNAL_ASSERT( B.columns() == (~lhs).columns() , "Invalid number of columns" );
1712 
1713  smpAddAssign( ~lhs, A * B );
1714  }
1716  //**********************************************************************************************
1717 
1718  //**Restructuring SMP addition assignment*******************************************************
1733  template< typename MT // Type of the target matrix
1734  , bool SO > // Storage order of the target matrix
1735  friend inline auto smpAddAssign( Matrix<MT,SO>& lhs, const SMatTDMatMultExpr& rhs )
1736  -> EnableIf_t< CanExploitSymmetry_v<MT,MT1,MT2> >
1737  {
1739 
1740  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1741  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1742 
1743  const ForwardFunctor fwd;
1744 
1745  smpAddAssign( ~lhs, fwd( rhs.lhs_ * trans( rhs.rhs_ ) ) );
1746  }
1748  //**********************************************************************************************
1749 
1750  //**SMP addition assignment to sparse matrices**************************************************
1751  // No special implementation for the SMP addition assignment to sparse matrices.
1752  //**********************************************************************************************
1753 
1754  //**SMP subtraction assignment to dense matrices************************************************
1770  template< typename MT // Type of the target dense matrix
1771  , bool SO > // Storage order of the target dense matrix
1772  friend inline auto smpSubAssign( DenseMatrix<MT,SO>& lhs, const SMatTDMatMultExpr& rhs )
1773  -> EnableIf_t< IsEvaluationRequired_v<MT,MT1,MT2> >
1774  {
1776 
1777  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1778  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1779 
1780  LT A( rhs.lhs_ ); // Evaluation of the left-hand side sparse matrix operand
1781  RT B( rhs.rhs_ ); // Evaluation of the right-hand side dense matrix operand
1782 
1783  BLAZE_INTERNAL_ASSERT( A.rows() == rhs.lhs_.rows() , "Invalid number of rows" );
1784  BLAZE_INTERNAL_ASSERT( A.columns() == rhs.lhs_.columns(), "Invalid number of columns" );
1785  BLAZE_INTERNAL_ASSERT( B.rows() == rhs.rhs_.rows() , "Invalid number of rows" );
1786  BLAZE_INTERNAL_ASSERT( B.columns() == rhs.rhs_.columns(), "Invalid number of columns" );
1787  BLAZE_INTERNAL_ASSERT( A.rows() == (~lhs).rows() , "Invalid number of rows" );
1788  BLAZE_INTERNAL_ASSERT( B.columns() == (~lhs).columns() , "Invalid number of columns" );
1789 
1790  smpSubAssign( ~lhs, A * B );
1791  }
1793  //**********************************************************************************************
1794 
1795  //**Restructuring SMP subtraction assignment****************************************************
1810  template< typename MT // Type of the target matrix
1811  , bool SO > // Storage order of the target matrix
1812  friend inline auto smpSubAssign( Matrix<MT,SO>& lhs, const SMatTDMatMultExpr& rhs )
1813  -> EnableIf_t< CanExploitSymmetry_v<MT,MT1,MT2> >
1814  {
1816 
1817  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1818  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1819 
1820  const ForwardFunctor fwd;
1821 
1822  smpSubAssign( ~lhs, fwd( rhs.lhs_ * trans( rhs.rhs_ ) ) );
1823  }
1825  //**********************************************************************************************
1826 
1827  //**SMP subtraction assignment to sparse matrices***********************************************
1828  // No special implementation for the SMP subtraction assignment to sparse matrices.
1829  //**********************************************************************************************
1830 
1831  //**SMP Schur product assignment to dense matrices**********************************************
1844  template< typename MT // Type of the target dense matrix
1845  , bool SO > // Storage order of the target dense matrix
1846  friend inline void smpSchurAssign( DenseMatrix<MT,SO>& lhs, const SMatTDMatMultExpr& rhs )
1847  {
1849 
1853 
1854  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1855  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1856 
1857  const ResultType tmp( rhs );
1858  smpSchurAssign( ~lhs, tmp );
1859  }
1861  //**********************************************************************************************
1862 
1863  //**SMP Schur product assignment to sparse matrices*********************************************
1864  // No special implementation for the SMP Schur product assignment to sparse matrices.
1865  //**********************************************************************************************
1866 
1867  //**SMP multiplication assignment to dense matrices*********************************************
1868  // No special implementation for the SMP multiplication assignment to dense matrices.
1869  //**********************************************************************************************
1870 
1871  //**SMP multiplication assignment to sparse matrices********************************************
1872  // No special implementation for the SMP multiplication assignment to sparse matrices.
1873  //**********************************************************************************************
1874 
1875  //**Compile time checks*************************************************************************
1884  //**********************************************************************************************
1885 };
1886 //*************************************************************************************************
1887 
1888 
1889 
1890 
1891 //=================================================================================================
1892 //
1893 // GLOBAL BINARY ARITHMETIC OPERATORS
1894 //
1895 //=================================================================================================
1896 
1897 //*************************************************************************************************
1910 template< typename MT1 // Type of the left-hand side dense matrix
1911  , typename MT2 // Type of the right-hand side sparse matrix
1912  , DisableIf_t< ( IsIdentity_v<MT1> &&
1913  IsSame_v< ElementType_t<MT1>, ElementType_t<MT2> > ) ||
1914  IsZero_v<MT1> >* = nullptr >
1915 inline const SMatTDMatMultExpr<MT1,MT2,false,false,false,false>
1916  smattdmatmult( const SparseMatrix<MT1,false>& lhs, const DenseMatrix<MT2,true>& rhs )
1917 {
1919 
1920  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == (~rhs).rows(), "Invalid matrix sizes" );
1921 
1922  return SMatTDMatMultExpr<MT1,MT2,false,false,false,false>( ~lhs, ~rhs );
1923 }
1925 //*************************************************************************************************
1926 
1927 
1928 //*************************************************************************************************
1942 template< typename MT1 // Type of the left-hand side sparse matrix
1943  , typename MT2 // Type of the right-hand side dense matrix
1944  , EnableIf_t< IsIdentity_v<MT1> &&
1945  IsSame_v< ElementType_t<MT1>, ElementType_t<MT2> > >* = nullptr >
1946 inline const MT2&
1947  smattdmatmult( const SparseMatrix<MT1,false>& lhs, const DenseMatrix<MT2,true>& rhs )
1948 {
1950 
1951  UNUSED_PARAMETER( lhs );
1952 
1953  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == (~rhs).rows(), "Invalid matrix sizes" );
1954 
1955  return (~rhs);
1956 }
1958 //*************************************************************************************************
1959 
1960 
1961 //*************************************************************************************************
1974 template< typename MT1 // Type of the left-hand side dense matrix
1975  , typename MT2 // Type of the right-hand side sparse matrix
1976  , EnableIf_t< IsZero_v<MT1> >* = nullptr >
1977 inline decltype(auto)
1978  smattdmatmult( const SparseMatrix<MT1,false>& lhs, const DenseMatrix<MT2,true>& rhs )
1979 {
1981 
1982  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == (~rhs).rows(), "Invalid matrix sizes" );
1983 
1984  using ReturnType = const MultTrait_t< ResultType_t<MT1>, ResultType_t<MT2> >;
1985 
1988 
1989  return ReturnType( (~lhs).rows(), (~rhs).columns() );
1990 }
1992 //*************************************************************************************************
1993 
1994 
1995 //*************************************************************************************************
2026 template< typename MT1 // Type of the left-hand side sparse matrix
2027  , typename MT2 > // Type of the right-hand side dense matrix
2028 inline decltype(auto)
2029  operator*( const SparseMatrix<MT1,false>& lhs, const DenseMatrix<MT2,true>& rhs )
2030 {
2032 
2033  if( (~lhs).columns() != (~rhs).rows() ) {
2034  BLAZE_THROW_INVALID_ARGUMENT( "Matrix sizes do not match" );
2035  }
2036 
2037  return smattdmatmult( ~lhs, ~rhs );
2038 }
2039 //*************************************************************************************************
2040 
2041 
2042 
2043 
2044 //=================================================================================================
2045 //
2046 // GLOBAL FUNCTIONS
2047 //
2048 //=================================================================================================
2049 
2050 //*************************************************************************************************
2076 template< typename MT1 // Type of the left-hand side dense matrix
2077  , typename MT2 // Type of the right-hand side dense matrix
2078  , bool SF // Symmetry flag
2079  , bool HF // Hermitian flag
2080  , bool LF // Lower flag
2081  , bool UF > // Upper flag
2082 inline decltype(auto) declsym( const SMatTDMatMultExpr<MT1,MT2,SF,HF,LF,UF>& dm )
2083 {
2085 
2086  if( !isSquare( dm ) ) {
2087  BLAZE_THROW_INVALID_ARGUMENT( "Invalid symmetric matrix specification" );
2088  }
2089 
2090  using ReturnType = const SMatTDMatMultExpr<MT1,MT2,true,HF,LF,UF>;
2091  return ReturnType( dm.leftOperand(), dm.rightOperand() );
2092 }
2094 //*************************************************************************************************
2095 
2096 
2097 //*************************************************************************************************
2123 template< typename MT1 // Type of the left-hand side dense matrix
2124  , typename MT2 // Type of the right-hand side dense matrix
2125  , bool SF // Symmetry flag
2126  , bool HF // Hermitian flag
2127  , bool LF // Lower flag
2128  , bool UF > // Upper flag
2129 inline decltype(auto) declherm( const SMatTDMatMultExpr<MT1,MT2,SF,HF,LF,UF>& dm )
2130 {
2132 
2133  if( !isSquare( dm ) ) {
2134  BLAZE_THROW_INVALID_ARGUMENT( "Invalid Hermitian matrix specification" );
2135  }
2136 
2137  using ReturnType = const SMatTDMatMultExpr<MT1,MT2,SF,true,LF,UF>;
2138  return ReturnType( dm.leftOperand(), dm.rightOperand() );
2139 }
2141 //*************************************************************************************************
2142 
2143 
2144 //*************************************************************************************************
2170 template< typename MT1 // Type of the left-hand side dense matrix
2171  , typename MT2 // Type of the right-hand side dense matrix
2172  , bool SF // Symmetry flag
2173  , bool HF // Hermitian flag
2174  , bool LF // Lower flag
2175  , bool UF > // Upper flag
2176 inline decltype(auto) decllow( const SMatTDMatMultExpr<MT1,MT2,SF,HF,LF,UF>& dm )
2177 {
2179 
2180  if( !isSquare( dm ) ) {
2181  BLAZE_THROW_INVALID_ARGUMENT( "Invalid lower matrix specification" );
2182  }
2183 
2184  using ReturnType = const SMatTDMatMultExpr<MT1,MT2,SF,HF,true,UF>;
2185  return ReturnType( dm.leftOperand(), dm.rightOperand() );
2186 }
2188 //*************************************************************************************************
2189 
2190 
2191 //*************************************************************************************************
2217 template< typename MT1 // Type of the left-hand side dense matrix
2218  , typename MT2 // Type of the right-hand side dense matrix
2219  , bool SF // Symmetry flag
2220  , bool HF // Hermitian flag
2221  , bool LF // Lower flag
2222  , bool UF > // Upper flag
2223 inline decltype(auto) declupp( const SMatTDMatMultExpr<MT1,MT2,SF,HF,LF,UF>& dm )
2224 {
2226 
2227  if( !isSquare( dm ) ) {
2228  BLAZE_THROW_INVALID_ARGUMENT( "Invalid upper matrix specification" );
2229  }
2230 
2231  using ReturnType = const SMatTDMatMultExpr<MT1,MT2,SF,HF,LF,true>;
2232  return ReturnType( dm.leftOperand(), dm.rightOperand() );
2233 }
2235 //*************************************************************************************************
2236 
2237 
2238 //*************************************************************************************************
2264 template< typename MT1 // Type of the left-hand side dense matrix
2265  , typename MT2 // Type of the right-hand side dense matrix
2266  , bool SF // Symmetry flag
2267  , bool HF // Hermitian flag
2268  , bool LF // Lower flag
2269  , bool UF > // Upper flag
2270 inline decltype(auto) decldiag( const SMatTDMatMultExpr<MT1,MT2,SF,HF,LF,UF>& dm )
2271 {
2273 
2274  if( !isSquare( dm ) ) {
2275  BLAZE_THROW_INVALID_ARGUMENT( "Invalid diagonal matrix specification" );
2276  }
2277 
2278  using ReturnType = const SMatTDMatMultExpr<MT1,MT2,SF,HF,true,true>;
2279  return ReturnType( dm.leftOperand(), dm.rightOperand() );
2280 }
2282 //*************************************************************************************************
2283 
2284 
2285 
2286 
2287 //=================================================================================================
2288 //
2289 // SIZE SPECIALIZATIONS
2290 //
2291 //=================================================================================================
2292 
2293 //*************************************************************************************************
2295 template< typename MT1, typename MT2, bool SF, bool HF, bool LF, bool UF >
2296 struct Size< SMatTDMatMultExpr<MT1,MT2,SF,HF,LF,UF>, 0UL >
2297  : public Size<MT1,0UL>
2298 {};
2299 
2300 template< typename MT1, typename MT2, bool SF, bool HF, bool LF, bool UF >
2301 struct Size< SMatTDMatMultExpr<MT1,MT2,SF,HF,LF,UF>, 1UL >
2302  : public Size<MT2,1UL>
2303 {};
2305 //*************************************************************************************************
2306 
2307 
2308 
2309 
2310 //=================================================================================================
2311 //
2312 // ISALIGNED SPECIALIZATIONS
2313 //
2314 //=================================================================================================
2315 
2316 //*************************************************************************************************
2318 template< typename MT1, typename MT2, bool SF, bool HF, bool LF, bool UF >
2319 struct IsAligned< SMatTDMatMultExpr<MT1,MT2,SF,HF,LF,UF> >
2320  : public IsAligned<MT2>
2321 {};
2323 //*************************************************************************************************
2324 
2325 } // namespace blaze
2326 
2327 #endif
decltype(auto) subvector(Vector< VT, TF > &, RSAs...)
Creating a view on a specific subvector of the given vector.
Definition: Subvector.h:329
#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.
decltype(auto) column(Matrix< MT, SO > &matrix, RCAs... args)
Creating a view on a specific column of the given matrix.
Definition: Column.h:133
Headerfile for the generic min algorithm.
Header file for the blaze::checked and blaze::unchecked instances.
Header file for the decldiag trait.
decltype(auto) decldiag(const DenseMatrix< MT, SO > &dm)
Declares the given dense matrix expression dm as diagonal.
Definition: DMatDeclDiagExpr.h:975
Header file for basic type definitions.
const ElementType ReturnType
Return type for expression template evaluations.
Definition: SMatTDMatMultExpr.h:236
typename If< Condition, T1, T2 >::Type If_t
Auxiliary alias declaration for the If class template.The If_t alias declaration provides a convenien...
Definition: If.h:109
RightOperand rightOperand() const noexcept
Returns the right-hand side transpose dense matrix operand.
Definition: SMatTDMatMultExpr.h:377
Header file for the declherm trait.
typename T::ResultType ResultType_t
Alias declaration for nested ResultType type definitions.The ResultType_t alias declaration provides ...
Definition: Aliases.h:390
static constexpr bool simdEnabled
Compilation switch for the expression template evaluation strategy.
Definition: SMatTDMatMultExpr.h:254
Header file for the serial shim.
ElementType_t< RT1 > ET1
Element type of the left-hand side dense matrix expression.
Definition: SMatTDMatMultExpr.h:132
Header file for the IsDiagonal type trait.
Base template for the DeclUppTrait class.
Definition: DeclUppTrait.h:134
#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 DeclUpp functor.
MT::Iterator begin(Matrix< MT, SO > &matrix, size_t i)
Returns an iterator to the first element of row/column i.
Definition: Matrix.h:372
If_t< IsExpression_v< MT1 >, const MT1, const MT1 &> LeftOperand
Composite type of the left-hand side sparse matrix expression.
Definition: SMatTDMatMultExpr.h:240
void reset(const DiagonalProxy< MT > &proxy)
Resetting the represented element to the default initial values.
Definition: DiagonalProxy.h:591
static constexpr bool smpAssignable
Compilation flag for SMP assignments.
Definition: CompressedMatrix.h:3113
size_t columns() const noexcept
Returns the current number of columns of the matrix.
Definition: SMatTDMatMultExpr.h:357
constexpr Unchecked unchecked
Global Unchecked instance.The blaze::unchecked instance is an optional token for the creation of view...
Definition: Check.h:138
Constraint on the data type.
Header file for the IsIdentity type trait.
OppositeType_t< ResultType > OppositeType
Result type with opposite storage order for expression template evaluations.
Definition: SMatTDMatMultExpr.h:233
decltype(auto) declupp(const DenseMatrix< MT, SO > &dm)
Declares the given dense matrix expression dm as upper.
Definition: DMatDeclUppExpr.h:1002
bool isAliased(const T *alias) const noexcept
Returns whether the expression is aliased with the given address alias.
Definition: SMatTDMatMultExpr.h:401
static constexpr bool evaluateRight
Compilation switch for the composite type of the right-hand side dense matrix expression.
Definition: SMatTDMatMultExpr.h:145
Header file for the Computation base class.
Header file for the MatMatMultExpr base class.
Header file for the reset shim.
Constraints on the storage order of matrix types.
CompositeType_t< MT1 > CT1
Composite type of the left-hand side sparse matrix expression.
Definition: SMatTDMatMultExpr.h:134
Header file for the RequiresEvaluation type trait.
typename If_t< HERM, DeclHermTrait< MultTrait_t< RT1, RT2 > >, If_t< SYM, DeclSymTrait< MultTrait_t< RT1, RT2 > >, If_t< LOW, If_t< UPP, DeclDiagTrait< MultTrait_t< RT1, RT2 > >, DeclLowTrait< MultTrait_t< RT1, RT2 > > >, If_t< UPP, DeclUppTrait< MultTrait_t< RT1, RT2 > >, MultTrait< RT1, RT2 > > > > >::Type ResultType
Result type for expression template evaluations.
Definition: SMatTDMatMultExpr.h:231
System settings for performance optimizations.
constexpr void UNUSED_PARAMETER(const Args &...)
Suppression of unused parameter warnings.
Definition: Unused.h:81
static constexpr bool evaluateLeft
Compilation switch for the composite type of the left-hand side sparse matrix expression.
Definition: SMatTDMatMultExpr.h:140
size_t rows() const noexcept
Returns the current number of rows of the matrix.
Definition: SMatTDMatMultExpr.h:347
constexpr size_t columns(const Matrix< MT, SO > &matrix) noexcept
Returns the current number of columns of the matrix.
Definition: Matrix.h:514
Base class for dense matrices.The DenseMatrix class is a base class for all dense matrix classes...
Definition: DenseMatrix.h:80
Base class for sparse matrices.The SparseMatrix class is a base class for all sparse matrix classes...
Definition: Forward.h:137
typename T::ElementType ElementType_t
Alias declaration for nested ElementType type definitions.The ElementType_t alias declaration provide...
Definition: Aliases.h:170
static constexpr bool UPP
Flag for upper matrices.
Definition: SMatTDMatMultExpr.h:152
Constraint on the data type.
Constraint on the data type.
bool canAlias(const T *alias) const noexcept
Returns whether the expression can alias with the given address alias.
Definition: SMatTDMatMultExpr.h:389
ElementType_t< ResultType > ElementType
Resulting element type.
Definition: SMatTDMatMultExpr.h:235
Headerfile for the generic max algorithm.
static constexpr bool HERM
Flag for Hermitian matrices.
Definition: SMatTDMatMultExpr.h:150
Header file for the DisableIf class template.
static constexpr bool SYM
Flag for symmetric matrices.
Definition: SMatTDMatMultExpr.h:149
Header file for the multiplication trait.
Header file for the IsStrictlyUpper type trait.
Header file for the IsSymmetric type trait.
Namespace of the Blaze C++ math library.
Definition: Blaze.h:58
Header file for the DeclLow functor.
Header file for the If class template.
#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
bool canSMPAssign() const noexcept
Returns whether the expression can be used in SMP assignments.
Definition: SMatTDMatMultExpr.h:421
#define BLAZE_CONSTRAINT_MUST_BE_ZERO_TYPE(T)
Constraint on the data type.In case the given data type T is not a zero vector or matrix type...
Definition: Zero.h:61
ResultType_t< MT2 > RT2
Result type of the right-hand side dense matrix expression.
Definition: SMatTDMatMultExpr.h:131
bool isAligned() const noexcept
Returns whether the operands of the expression are properly aligned in memory.
Definition: SMatTDMatMultExpr.h:411
Generic wrapper for the decllow() function.
Definition: DeclLow.h:58
CompositeType_t< MT2 > CT2
Composite type of the right-hand side dense matrix expression.
Definition: SMatTDMatMultExpr.h:135
decltype(auto) min(const DenseMatrix< MT1, SO1 > &lhs, const DenseMatrix< MT2, SO2 > &rhs)
Computes the componentwise minimum of the dense matrices lhs and rhs.
Definition: DMatDMatMapExpr.h:1147
Header file for the decllow trait.
#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 DenseMatrix base class.
If_t< evaluateLeft, const RT1, CT1 > LT
Type for the assignment of the left-hand side sparse matrix operand.
Definition: SMatTDMatMultExpr.h:246
decltype(auto) decllow(const DenseMatrix< MT, SO > &dm)
Declares the given dense matrix expression dm as lower.
Definition: DMatDeclLowExpr.h:1002
Header file for the IsLower type trait.
LeftOperand lhs_
Left-hand side sparse matrix of the multiplication expression.
Definition: SMatTDMatMultExpr.h:428
Header file for the IsAligned type trait.
Expression object for sparse matrix-transpose dense matrix multiplications.The SMatTDMatMultExpr clas...
Definition: Forward.h:129
const ResultType CompositeType
Data type for composite expression templates.
Definition: SMatTDMatMultExpr.h:237
Generic wrapper for the null function.
Definition: Noop.h:59
Header file for the IsTriangular type trait.
Base template for the DeclSymTrait class.
Definition: DeclSymTrait.h:134
Constraints on the storage order of matrix types.
Header file for the exception macros of the math module.
decltype(auto) max(const DenseMatrix< MT1, SO1 > &lhs, const DenseMatrix< MT2, SO2 > &rhs)
Computes the componentwise maximum of the dense matrices lhs and rhs.
Definition: DMatDMatMapExpr.h:1179
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:438
Header file for the DeclDiag functor.
Constraint on the data type.
Header file for all forward declarations for expression class templates.
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:103
typename T::OppositeType OppositeType_t
Alias declaration for nested OppositeType type definitions.The OppositeType_t alias declaration provi...
Definition: Aliases.h:270
Header file for the conjugate shim.
Header file for the declupp trait.
#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
static constexpr bool LOW
Flag for lower matrices.
Definition: SMatTDMatMultExpr.h:151
typename T::TransposeType TransposeType_t
Alias declaration for nested TransposeType type definitions.The TransposeType_t alias declaration pro...
Definition: Aliases.h:470
Header file for run time assertion macros.
Base template for the DeclHermTrait class.
Definition: DeclHermTrait.h:134
typename T::CompositeType CompositeType_t
Alias declaration for nested CompositeType type definitions.The CompositeType_t alias declaration pro...
Definition: Aliases.h:90
Base template for the MultTrait class.
Definition: MultTrait.h:146
auto smpAddAssign(Matrix< MT1, SO1 > &lhs, const Matrix< MT2, SO2 > &rhs) -> EnableIf_t< IsDenseMatrix_v< MT1 > >
Default implementation of the SMP addition assignment of a matrix to a dense matrix.
Definition: DenseMatrix.h:131
static constexpr bool smpAssignable
Compilation switch for the expression template assignment strategy.
Definition: SMatTDMatMultExpr.h:257
decltype(auto) row(Matrix< MT, SO > &, RRAs...)
Creating a view on a specific row of the given matrix.
Definition: Row.h:133
TransposeType_t< ResultType > TransposeType
Transpose type for expression template evaluations.
Definition: SMatTDMatMultExpr.h:234
Header file for the IsZero type trait.
Header file for the declsym trait.
#define BLAZE_FUNCTION_TRACE
Function trace macro.This macro can be used to reliably trace function calls. In case function tracin...
Definition: FunctionTrace.h:94
decltype(auto) declsym(const DenseMatrix< MT, SO > &dm)
Declares the given dense matrix expression dm as symmetric.
Definition: DMatDeclSymExpr.h:1002
auto smpAssign(Matrix< MT1, SO1 > &lhs, const Matrix< MT2, SO2 > &rhs) -> EnableIf_t< IsDenseMatrix_v< MT1 > >
Default implementation of the SMP assignment of a matrix to a dense matrix.
Definition: DenseMatrix.h:100
Constraints on the storage order of matrix types.
Generic wrapper for the declherm() function.
Definition: DeclHerm.h:58
decltype(auto) serial(const DenseMatrix< MT, SO > &dm)
Forces the serial evaluation of the given dense matrix expression dm.
Definition: DMatSerialExpr.h:808
SMatTDMatMultExpr(const MT1 &lhs, const MT2 &rhs) noexcept
Constructor for the SMatTDMatMultExpr class.
Definition: SMatTDMatMultExpr.h:267
Header file for the Noop functor.
#define BLAZE_CONSTRAINT_MUST_NOT_REQUIRE_EVALUATION(T)
Constraint on the data type.In case the given data type T requires an intermediate evaluation within ...
Definition: RequiresEvaluation.h:81
auto smpSchurAssign(Matrix< MT1, SO1 > &lhs, const Matrix< MT2, SO2 > &rhs) -> EnableIf_t< IsDenseMatrix_v< MT1 > >
Default implementation of the SMP Schur product assignment of a matrix to dense matrix.
Definition: DenseMatrix.h:194
constexpr size_t rows(const Matrix< MT, SO > &matrix) noexcept
Returns the current number of rows of the matrix.
Definition: Matrix.h:498
#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
Generic wrapper for the declupp() function.
Definition: DeclUpp.h:58
const Type & ReturnType
Return type for expression template evaluations.
Definition: CompressedMatrix.h:3081
decltype(auto) trans(const DenseMatrix< MT, SO > &dm)
Calculation of the transpose of the given dense matrix.
Definition: DMatTransExpr.h:765
LeftOperand leftOperand() const noexcept
Returns the left-hand side sparse matrix operand.
Definition: SMatTDMatMultExpr.h:367
Base template for the DeclLowTrait class.
Definition: DeclLowTrait.h:134
decltype(auto) declherm(const DenseMatrix< MT, SO > &dm)
Declares the given dense matrix expression dm as Hermitian.
Definition: DMatDeclHermExpr.h:1002
Header file for the IsComputation type trait class.
RightOperand rhs_
Right-hand side dense matrix of the multiplication expression.
Definition: SMatTDMatMultExpr.h:429
Header file for the IsBuiltin type trait.
auto smpSubAssign(Matrix< MT1, SO1 > &lhs, const Matrix< MT2, SO2 > &rhs) -> EnableIf_t< IsDenseMatrix_v< MT1 > >
Default implementation of the SMP subtraction assignment of a matrix to dense matrix.
Definition: DenseMatrix.h:162
ReturnType operator()(size_t i, size_t j) const
2D-access to the matrix elements.
Definition: SMatTDMatMultExpr.h:282
ResultType_t< MT1 > RT1
Result type of the left-hand side sparse matrix expression.
Definition: SMatTDMatMultExpr.h:130
Header file for the IntegralConstant class template.
Generic wrapper for the decldiag() function.
Definition: DeclDiag.h:58
Header file for the DeclHerm functor.
If_t< evaluateRight, const RT2, CT2 > RT
Type for the assignment of the right-hand side dense matrix operand.
Definition: SMatTDMatMultExpr.h:249
If_t< IsExpression_v< MT2 >, const MT2, const MT2 &> RightOperand
Composite type of the right-hand side dense matrix expression.
Definition: SMatTDMatMultExpr.h:243
ElementType_t< RT2 > ET2
Element type of the right-hand side sparse matrix expression.
Definition: SMatTDMatMultExpr.h:133
Header file for the IsUpper type trait.
typename DisableIf< Condition, T >::Type DisableIf_t
Auxiliary type for the DisableIf class template.The DisableIf_t alias declaration provides a convenie...
Definition: DisableIf.h:138
decltype(auto) conj(const DenseMatrix< MT, SO > &dm)
Returns a matrix containing the complex conjugate of each single element of dm.
Definition: DMatMapExpr.h:1326
Constraint on the data type.
Generic wrapper for the declsym() function.
Definition: DeclSym.h:58
Base template for the DeclDiagTrait class.
Definition: DeclDiagTrait.h:134
bool isSquare(const Matrix< MT, SO > &matrix) noexcept
Checks if the given matrix is a square matrix.
Definition: Matrix.h:951
Header file for the IsResizable type trait.
Header file for the Size type trait.
#define BLAZE_CONSTRAINT_MUST_NOT_BE_ZERO_TYPE(T)
Constraint on the data type.In case the given data type T is a zero vector or matrix type...
Definition: Zero.h:81
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
Header file for the DeclSym functor.
#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
Header file for the TrueType type/value trait base class.
Header file for the IsExpression type trait class.
Header file for the function trace functionality.
ReturnType at(size_t i, size_t j) const
Checked access to the matrix elements.
Definition: SMatTDMatMultExpr.h:331