DMatTSMatMultExpr.h
Go to the documentation of this file.
1 //=================================================================================================
33 //=================================================================================================
34 
35 #ifndef _BLAZE_MATH_EXPRESSIONS_DMATTSMATMULTEXPR_H_
36 #define _BLAZE_MATH_EXPRESSIONS_DMATTSMATMULTEXPR_H_
37 
38 
39 //*************************************************************************************************
40 // Includes
41 //*************************************************************************************************
42 
43 #include <blaze/math/Aliases.h>
53 #include <blaze/math/Exception.h>
65 #include <blaze/math/shims/Reset.h>
88 #include <blaze/math/views/Check.h>
93 #include <blaze/util/Assert.h>
94 #include <blaze/util/DisableIf.h>
95 #include <blaze/util/EnableIf.h>
98 #include <blaze/util/mpl/If.h>
99 #include <blaze/util/TrueType.h>
100 #include <blaze/util/Types.h>
103 #include <blaze/util/Unused.h>
104 
105 
106 namespace blaze {
107 
108 //=================================================================================================
109 //
110 // CLASS DMATTSMATMULTEXPR
111 //
112 //=================================================================================================
113 
114 //*************************************************************************************************
121 template< typename MT1 // Type of the left-hand side dense matrix
122  , typename MT2 // Type of the right-hand side sparse matrix
123  , bool SF // Symmetry flag
124  , bool HF // Hermitian flag
125  , bool LF // Lower flag
126  , bool UF > // Upper flag
128  : public MatMatMultExpr< DenseMatrix< DMatTSMatMultExpr<MT1,MT2,SF,HF,LF,UF>, false > >
129  , private Computation
130 {
131  private:
132  //**Type definitions****************************************************************************
139  //**********************************************************************************************
140 
141  //**********************************************************************************************
143  static constexpr bool evaluateLeft = ( IsComputation_v<MT1> || RequiresEvaluation_v<MT1> );
144  //**********************************************************************************************
145 
146  //**********************************************************************************************
148  static constexpr bool evaluateRight = ( IsComputation_v<MT2> || RequiresEvaluation_v<MT2> );
149  //**********************************************************************************************
150 
151  //**********************************************************************************************
152  static constexpr bool SYM = ( SF && !( HF || LF || UF ) );
153  static constexpr bool HERM = ( HF && !( LF || UF ) );
154  static constexpr bool LOW = ( LF || ( ( SF || HF ) && UF ) );
155  static constexpr bool UPP = ( UF || ( ( SF || HF ) && LF ) );
156  //**********************************************************************************************
157 
158  //**********************************************************************************************
160 
165  template< typename T1, typename T2, typename T3 >
166  static constexpr bool CanExploitSymmetry_v = IsSymmetric_v<T2>;
168  //**********************************************************************************************
169 
170  //**********************************************************************************************
172 
176  template< typename T1, typename T2, typename T3 >
177  static constexpr bool IsEvaluationRequired_v =
178  ( ( evaluateLeft || evaluateRight ) && !CanExploitSymmetry_v<T1,T2,T3> );
180  //**********************************************************************************************
181 
182  //**********************************************************************************************
184 
187  template< typename T1, typename T2, typename T3 >
188  static constexpr bool UseOptimizedKernel_v =
189  ( useOptimizedKernels &&
190  !IsDiagonal_v<T2> &&
191  !IsResizable_v< ElementType_t<T1> > &&
192  !IsResizable_v<ET2> );
194  //**********************************************************************************************
195 
196  //**********************************************************************************************
198 
201  using ForwardFunctor = If_t< HERM
202  , DeclHerm
203  , If_t< SYM
204  , DeclSym
205  , If_t< LOW
206  , If_t< UPP
207  , DeclDiag
208  , DeclLow >
209  , If_t< UPP
210  , DeclUpp
211  , Noop > > > >;
213  //**********************************************************************************************
214 
215  public:
216  //**Type definitions****************************************************************************
219 
222 
224  using ResultType = typename If_t< HERM
226  , If_t< SYM
228  , If_t< LOW
229  , If_t< UPP
232  , If_t< UPP
234  , MultTrait<RT1,RT2> > > > >::Type;
235 
239  using ReturnType = const ElementType;
240  using CompositeType = const ResultType;
241 
243  using LeftOperand = If_t< IsExpression_v<MT1>, const MT1, const MT1& >;
244 
246  using RightOperand = If_t< IsExpression_v<MT2>, const MT2, const MT2& >;
247 
250 
253  //**********************************************************************************************
254 
255  //**Compilation flags***************************************************************************
257  static constexpr bool simdEnabled = false;
258 
260  static constexpr bool smpAssignable =
262  //**********************************************************************************************
263 
264  //**Constructor*********************************************************************************
270  explicit inline DMatTSMatMultExpr( const MT1& lhs, const MT2& rhs ) noexcept
271  : lhs_( lhs ) // Left-hand side dense matrix of the multiplication expression
272  , rhs_( rhs ) // Right-hand side sparse matrix of the multiplication expression
273  {
274  BLAZE_INTERNAL_ASSERT( lhs.columns() == rhs.rows(), "Invalid matrix sizes" );
275  }
276  //**********************************************************************************************
277 
278  //**Access operator*****************************************************************************
285  inline ReturnType operator()( size_t i, size_t j ) const {
286  BLAZE_INTERNAL_ASSERT( i < lhs_.rows() , "Invalid row access index" );
287  BLAZE_INTERNAL_ASSERT( j < rhs_.columns(), "Invalid column access index" );
288 
289  if( IsDiagonal_v<MT1> ) {
290  return lhs_(i,i) * rhs_(i,j);
291  }
292  else if( IsDiagonal_v<MT2> ) {
293  return lhs_(i,j) * rhs_(j,j);
294  }
295  else if( IsTriangular_v<MT1> || IsTriangular_v<MT2> ) {
296  const size_t begin( ( IsUpper_v<MT1> )
297  ?( ( IsLower_v<MT2> )
298  ?( max( ( IsStrictlyUpper_v<MT1> ? i+1UL : i )
299  , ( IsStrictlyLower_v<MT2> ? j+1UL : j ) ) )
300  :( IsStrictlyUpper_v<MT1> ? i+1UL : i ) )
301  :( ( IsLower_v<MT2> )
302  ?( IsStrictlyLower_v<MT2> ? j+1UL : j )
303  :( 0UL ) ) );
304  const size_t end( ( IsLower_v<MT1> )
305  ?( ( IsUpper_v<MT2> )
306  ?( min( ( IsStrictlyLower_v<MT1> ? i : i+1UL )
307  , ( IsStrictlyUpper_v<MT2> ? j : j+1UL ) ) )
308  :( IsStrictlyLower_v<MT1> ? i : i+1UL ) )
309  :( ( IsUpper_v<MT2> )
310  ?( IsStrictlyUpper_v<MT2> ? j : j+1UL )
311  :( lhs_.columns() ) ) );
312 
313  if( begin >= end ) return ElementType();
314 
315  const size_t n( end - begin );
316 
317  return subvector( row( lhs_, i, unchecked ), begin, n, unchecked ) *
318  subvector( column( rhs_, j, unchecked ), begin, n, unchecked );
319  }
320  else {
321  return row( lhs_, i, unchecked ) * column( rhs_, j, unchecked );
322  }
323  }
324  //**********************************************************************************************
325 
326  //**At function*********************************************************************************
334  inline ReturnType at( size_t i, size_t j ) const {
335  if( i >= lhs_.rows() ) {
336  BLAZE_THROW_OUT_OF_RANGE( "Invalid row access index" );
337  }
338  if( j >= rhs_.columns() ) {
339  BLAZE_THROW_OUT_OF_RANGE( "Invalid column access index" );
340  }
341  return (*this)(i,j);
342  }
343  //**********************************************************************************************
344 
345  //**Rows function*******************************************************************************
350  inline size_t rows() const noexcept {
351  return lhs_.rows();
352  }
353  //**********************************************************************************************
354 
355  //**Columns function****************************************************************************
360  inline size_t columns() const noexcept {
361  return rhs_.columns();
362  }
363  //**********************************************************************************************
364 
365  //**Left operand access*************************************************************************
370  inline LeftOperand leftOperand() const noexcept {
371  return lhs_;
372  }
373  //**********************************************************************************************
374 
375  //**Right operand access************************************************************************
380  inline RightOperand rightOperand() const noexcept {
381  return rhs_;
382  }
383  //**********************************************************************************************
384 
385  //**********************************************************************************************
391  template< typename T >
392  inline bool canAlias( const T* alias ) const noexcept {
393  return ( lhs_.isAliased( alias ) || rhs_.isAliased( alias ) );
394  }
395  //**********************************************************************************************
396 
397  //**********************************************************************************************
403  template< typename T >
404  inline bool isAliased( const T* alias ) const noexcept {
405  return ( lhs_.isAliased( alias ) || rhs_.isAliased( alias ) );
406  }
407  //**********************************************************************************************
408 
409  //**********************************************************************************************
414  inline bool isAligned() const noexcept {
415  return lhs_.isAligned();
416  }
417  //**********************************************************************************************
418 
419  //**********************************************************************************************
424  inline bool canSMPAssign() const noexcept {
425  return ( rows() * columns() >= SMP_DMATTSMATMULT_THRESHOLD ) && !IsDiagonal_v<MT1>;
426  }
427  //**********************************************************************************************
428 
429  private:
430  //**Member variables****************************************************************************
433  //**********************************************************************************************
434 
435  //**Assignment to dense matrices****************************************************************
448  template< typename MT // Type of the target dense matrix
449  , bool SO > // Storage order of the target dense matrix
450  friend inline auto assign( DenseMatrix<MT,SO>& lhs, const DMatTSMatMultExpr& rhs )
452  {
454 
455  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
456  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
457 
458  LT A( serial( rhs.lhs_ ) ); // Evaluation of the left-hand side dense matrix operand
459  RT B( serial( rhs.rhs_ ) ); // Evaluation of the right-hand side sparse matrix operand
460 
461  BLAZE_INTERNAL_ASSERT( A.rows() == (~lhs).rows() , "Invalid number of rows" );
462  BLAZE_INTERNAL_ASSERT( A.columns() == B.rows() , "Invalid matrix sizes" );
463  BLAZE_INTERNAL_ASSERT( B.columns() == (~lhs).columns(), "Invalid number of columns" );
464 
465  DMatTSMatMultExpr::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  if( LOW && UPP ) {
496  reset( C );
497  }
498 
499  {
500  size_t i( 0UL );
501 
502  for( ; (i+4UL) <= M; i+=4UL ) {
503  for( size_t j=( SYM || HERM || UPP ? i : 0UL ); j<( LOW ? i+4UL : N ); ++j )
504  {
505  auto element( ( IsUpper_v<MT4> )
506  ?( IsStrictlyUpper_v<MT4> ? B.upperBound(i,j) : B.lowerBound(i,j) )
507  :( B.begin(j) ) );
508  const auto end( ( IsLower_v<MT4> )
509  ?( IsStrictlyLower_v<MT4> ? B.lowerBound(i+4UL,j) : B.upperBound(i+4UL,j) )
510  :( B.end(j) ) );
511 
512  if( element == end ) {
513  reset( C(i ,j) );
514  reset( C(i+1UL,j) );
515  reset( C(i+2UL,j) );
516  reset( C(i+3UL,j) );
517  continue;
518  }
519 
520  C(i ,j) = A(i ,element->index()) * element->value();
521  C(i+1UL,j) = A(i+1UL,element->index()) * element->value();
522  C(i+2UL,j) = A(i+2UL,element->index()) * element->value();
523  C(i+3UL,j) = A(i+3UL,element->index()) * element->value();
524  ++element;
525  for( ; element!=end; ++element ) {
526  C(i ,j) += A(i ,element->index()) * element->value();
527  C(i+1UL,j) += A(i+1UL,element->index()) * element->value();
528  C(i+2UL,j) += A(i+2UL,element->index()) * element->value();
529  C(i+3UL,j) += A(i+3UL,element->index()) * element->value();
530  }
531  }
532  }
533 
534  for( ; (i+2UL) <= M; i+=2UL ) {
535  for( size_t j=( SYM || HERM || UPP ? i : 0UL ); j<( LOW ? i+2UL : N ); ++j )
536  {
537  auto element( ( IsUpper_v<MT4> )
538  ?( IsStrictlyUpper_v<MT4> ? B.upperBound(i,j) : B.lowerBound(i,j) )
539  :( B.begin(j) ) );
540  const auto end( ( IsLower_v<MT4> )
541  ?( IsStrictlyLower_v<MT4> ? B.lowerBound(i+2UL,j) : B.upperBound(i+2UL,j) )
542  :( B.end(j) ) );
543 
544  if( element == end ) {
545  reset( C(i ,j) );
546  reset( C(i+1UL,j) );
547  continue;
548  }
549 
550  C(i ,j) = A(i ,element->index()) * element->value();
551  C(i+1UL,j) = A(i+1UL,element->index()) * element->value();
552  ++element;
553  for( ; element!=end; ++element ) {
554  C(i ,j) += A(i ,element->index()) * element->value();
555  C(i+1UL,j) += A(i+1UL,element->index()) * element->value();
556  }
557  }
558  }
559 
560  for( ; i<M; ++i ) {
561  for( size_t j=( SYM || HERM || UPP ? i : 0UL ); j<( LOW ? i+1UL : N ); ++j )
562  {
563  auto element( ( IsUpper_v<MT4> )
564  ?( IsStrictlyUpper_v<MT4> ? B.upperBound(i,j) : B.lowerBound(i,j) )
565  :( B.begin(j) ) );
566  const auto end( ( IsLower_v<MT4> )
567  ?( IsStrictlyLower_v<MT4> ? B.lowerBound(i,j) : B.upperBound(i,j) )
568  :( B.end(j) ) );
569 
570  if( element == end ) {
571  reset( C(i,j) );
572  continue;
573  }
574 
575  C(i,j) = A(i,element->index()) * element->value();
576  ++element;
577  for( ; element!=end; ++element )
578  C(i,j) += A(i,element->index()) * element->value();
579  }
580  }
581  }
582 
583  if( SYM || HERM ) {
584  for( size_t i=1UL; i<M; ++i ) {
585  for( size_t j=0UL; j<i; ++j ) {
586  C(i,j) = HERM ? conj( C(j,i) ) : C(j,i);
587  }
588  }
589  }
590  else if( LOW && !UPP ) {
591  for( size_t j=1UL; j<N; ++j ) {
592  for( size_t i=0UL; i<j; ++i ) {
593  reset( C(i,j) );
594  }
595  }
596  }
597  else if( !LOW && UPP ) {
598  for( size_t i=1UL; i<M; ++i ) {
599  for( size_t j=0UL; j<i; ++j ) {
600  reset( C(i,j) );
601  }
602  }
603  }
604  }
606  //**********************************************************************************************
607 
608  //**Optimized assignment to dense matrices******************************************************
622  template< typename MT3 // Type of the left-hand side target matrix
623  , typename MT4 // Type of the left-hand side matrix operand
624  , typename MT5 > // Type of the right-hand side matrix operand
625  static inline auto selectAssignKernel( MT3& C, const MT4& A, const MT5& B )
626  -> EnableIf_t< UseOptimizedKernel_v<MT3,MT4,MT5> >
627  {
628  const size_t M( A.rows() );
629  const size_t N( B.columns() );
630 
631  BLAZE_INTERNAL_ASSERT( !( SYM || HERM || LOW || UPP ) || M == N, "Broken invariant detected" );
632 
633  reset( C );
634 
635  {
636  size_t i( 0UL );
637 
638  for( ; (i+4UL) <= M; i+=4UL ) {
639  for( size_t j=( SYM || HERM || UPP ? i : 0UL ); j<( LOW ? i+4UL : N ); ++j )
640  {
641  auto element( ( IsUpper_v<MT4> )
642  ?( IsStrictlyUpper_v<MT4> ? B.upperBound(i,j) : B.lowerBound(i,j) )
643  :( B.begin(j) ) );
644  const auto end( ( IsLower_v<MT4> )
645  ?( IsStrictlyLower_v<MT4> ? B.lowerBound(i+4UL,j) : B.upperBound(i+4UL,j) )
646  :( B.end(j) ) );
647 
648  const size_t nonzeros( end - element );
649  const size_t kpos( nonzeros & size_t(-4) );
650  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
651 
652  for( size_t k=0UL; k<kpos; k+=4UL )
653  {
654  const size_t j1( element->index() );
655  const ET2 v1( element->value() );
656  ++element;
657  const size_t j2( element->index() );
658  const ET2 v2( element->value() );
659  ++element;
660  const size_t j3( element->index() );
661  const ET2 v3( element->value() );
662  ++element;
663  const size_t j4( element->index() );
664  const ET2 v4( element->value() );
665  ++element;
666 
667  BLAZE_INTERNAL_ASSERT( j1 < j2 && j2 < j3 && j3 < j4, "Invalid sparse matrix index detected" );
668 
669  C(i ,j) += A(i ,j1) * v1 + A(i ,j2) * v2 + A(i ,j3) * v3 + A(i ,j4) * v4;
670  C(i+1UL,j) += A(i+1UL,j1) * v1 + A(i+1UL,j2) * v2 + A(i+1UL,j3) * v3 + A(i+1UL,j4) * v4;
671  C(i+2UL,j) += A(i+2UL,j1) * v1 + A(i+2UL,j2) * v2 + A(i+2UL,j3) * v3 + A(i+2UL,j4) * v4;
672  C(i+3UL,j) += A(i+3UL,j1) * v1 + A(i+3UL,j2) * v2 + A(i+3UL,j3) * v3 + A(i+3UL,j4) * v4;
673  }
674 
675  for( ; element!=end; ++element )
676  {
677  const size_t j1( element->index() );
678  const ET2 v1( element->value() );
679 
680  C(i ,j) += A(i ,j1) * v1;
681  C(i+1UL,j) += A(i+1UL,j1) * v1;
682  C(i+2UL,j) += A(i+2UL,j1) * v1;
683  C(i+3UL,j) += A(i+3UL,j1) * v1;
684  }
685  }
686  }
687 
688  for( ; (i+2UL) <= M; i+=2UL ) {
689  for( size_t j=( SYM || HERM || UPP ? i : 0UL ); j<( LOW ? i+2UL : N ); ++j )
690  {
691  auto element( ( IsUpper_v<MT4> )
692  ?( IsStrictlyUpper_v<MT4> ? B.upperBound(i,j) : B.lowerBound(i,j) )
693  :( B.begin(j) ) );
694  const auto end( ( IsLower_v<MT4> )
695  ?( IsStrictlyLower_v<MT4> ? B.lowerBound(i+2UL,j) : B.upperBound(i+2UL,j) )
696  :( B.end(j) ) );
697 
698  const size_t nonzeros( end - element );
699  const size_t kpos( nonzeros & size_t(-4) );
700  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
701 
702  for( size_t k=0UL; k<kpos; k+=4UL )
703  {
704  const size_t j1( element->index() );
705  const ET2 v1( element->value() );
706  ++element;
707  const size_t j2( element->index() );
708  const ET2 v2( element->value() );
709  ++element;
710  const size_t j3( element->index() );
711  const ET2 v3( element->value() );
712  ++element;
713  const size_t j4( element->index() );
714  const ET2 v4( element->value() );
715  ++element;
716 
717  BLAZE_INTERNAL_ASSERT( j1 < j2 && j2 < j3 && j3 < j4, "Invalid sparse matrix index detected" );
718 
719  C(i ,j) += A(i ,j1) * v1 + A(i ,j2) * v2 + A(i ,j3) * v3 + A(i ,j4) * v4;
720  C(i+1UL,j) += A(i+1UL,j1) * v1 + A(i+1UL,j2) * v2 + A(i+1UL,j3) * v3 + A(i+1UL,j4) * v4;
721  }
722 
723  for( ; element!=end; ++element )
724  {
725  const size_t j1( element->index() );
726  const ET2 v1( element->value() );
727 
728  C(i ,j) += A(i ,j1) * v1;
729  C(i+1UL,j) += A(i+1UL,j1) * v1;
730  }
731  }
732  }
733 
734  for( ; i<M; ++i ) {
735  for( size_t j=( SYM || HERM || UPP ? i : 0UL ); j<( LOW ? i+1UL : N ); ++j )
736  {
737  auto element( ( IsUpper_v<MT4> )
738  ?( IsStrictlyUpper_v<MT4> ? B.upperBound(i,j) : B.lowerBound(i,j) )
739  :( B.begin(j) ) );
740  const auto end( ( IsLower_v<MT4> )
741  ?( IsStrictlyLower_v<MT4> ? B.lowerBound(i,j) : B.upperBound(i,j) )
742  :( B.end(j) ) );
743 
744  const size_t nonzeros( end - element );
745  const size_t kpos( nonzeros & size_t(-4) );
746  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
747 
748  for( size_t k=0UL; k<kpos; k+=4UL )
749  {
750  const size_t j1( element->index() );
751  const ET2 v1( element->value() );
752  ++element;
753  const size_t j2( element->index() );
754  const ET2 v2( element->value() );
755  ++element;
756  const size_t j3( element->index() );
757  const ET2 v3( element->value() );
758  ++element;
759  const size_t j4( element->index() );
760  const ET2 v4( element->value() );
761  ++element;
762 
763  BLAZE_INTERNAL_ASSERT( j1 < j2 && j2 < j3 && j3 < j4, "Invalid sparse matrix index detected" );
764 
765  C(i,j) += A(i,j1) * v1 + A(i,j2) * v2 + A(i,j3) * v3 + A(i,j4) * v4;
766  }
767 
768  for( ; element!=end; ++element )
769  {
770  const size_t j1( element->index() );
771  const ET2 v1( element->value() );
772 
773  C(i,j) += A(i,j1) * v1;
774  }
775  }
776  }
777  }
778 
779  if( SYM || HERM ) {
780  for( size_t i=1UL; i<M; ++i ) {
781  for( size_t j=0UL; j<i; ++j ) {
782  C(i,j) = HERM ? conj( C(j,i) ) : C(j,i);
783  }
784  }
785  }
786  }
788  //**********************************************************************************************
789 
790  //**Assignment to sparse matrices***************************************************************
803  template< typename MT // Type of the target sparse matrix
804  , bool SO > // Storage order of the target sparse matrix
805  friend inline auto assign( SparseMatrix<MT,SO>& lhs, const DMatTSMatMultExpr& rhs )
806  -> DisableIf_t< CanExploitSymmetry_v<MT,MT1,MT2> >
807  {
809 
810  using TmpType = If_t< SO, OppositeType, ResultType >;
811 
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  const ForwardFunctor fwd;
823 
824  const TmpType tmp( serial( rhs ) );
825  assign( ~lhs, fwd( tmp ) );
826  }
828  //**********************************************************************************************
829 
830  //**Restructuring assignment********************************************************************
845  template< typename MT // Type of the target matrix
846  , bool SO > // Storage order of the target matrix
847  friend inline auto assign( Matrix<MT,SO>& lhs, const DMatTSMatMultExpr& rhs )
848  -> EnableIf_t< CanExploitSymmetry_v<MT,MT1,MT2> >
849  {
851 
852  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
853  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
854 
855  const ForwardFunctor fwd;
856 
857  assign( ~lhs, fwd( trans( rhs.lhs_ ) * rhs.rhs_ ) );
858  }
860  //**********************************************************************************************
861 
862  //**Addition assignment to dense matrices*******************************************************
875  template< typename MT // Type of the target dense matrix
876  , bool SO > // Storage order of the target dense matrix
877  friend inline auto addAssign( DenseMatrix<MT,SO>& lhs, const DMatTSMatMultExpr& rhs )
878  -> DisableIf_t< CanExploitSymmetry_v<MT,MT1,MT2> >
879  {
881 
882  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
883  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
884 
885  LT A( serial( rhs.lhs_ ) ); // Evaluation of the left-hand side dense matrix operand
886  RT B( serial( rhs.rhs_ ) ); // Evaluation of the right-hand side sparse matrix operand
887 
888  BLAZE_INTERNAL_ASSERT( A.rows() == (~lhs).rows() , "Invalid number of rows" );
889  BLAZE_INTERNAL_ASSERT( A.columns() == B.rows() , "Invalid matrix sizes" );
890  BLAZE_INTERNAL_ASSERT( B.columns() == (~lhs).columns(), "Invalid number of columns" );
891 
892  DMatTSMatMultExpr::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 i( 0UL );
924 
925  for( ; (i+4UL) <= M; i+=4UL ) {
926  for( size_t j=( UPP ? i : 0UL ); j<( LOW ? i+4UL : N ); ++j )
927  {
928  auto element( ( IsUpper_v<MT4> )
929  ?( IsStrictlyUpper_v<MT4> ? B.upperBound(i,j) : B.lowerBound(i,j) )
930  :( B.begin(j) ) );
931  const auto end( ( IsLower_v<MT4> )
932  ?( IsStrictlyLower_v<MT4> ? B.lowerBound(i+4UL,j) : B.upperBound(i+4UL,j) )
933  :( B.end(j) ) );
934 
935  for( ; element!=end; ++element ) {
936  C(i ,j) += A(i ,element->index()) * element->value();
937  C(i+1UL,j) += A(i+1UL,element->index()) * element->value();
938  C(i+2UL,j) += A(i+2UL,element->index()) * element->value();
939  C(i+3UL,j) += A(i+3UL,element->index()) * element->value();
940  }
941  }
942  }
943 
944  for( ; (i+2UL) <= M; i+=2UL ) {
945  for( size_t j=( UPP ? i : 0UL ); j<( LOW ? i+2UL : N ); ++j )
946  {
947  auto element( ( IsUpper_v<MT4> )
948  ?( IsStrictlyUpper_v<MT4> ? B.upperBound(i,j) : B.lowerBound(i,j) )
949  :( B.begin(j) ) );
950  const auto end( ( IsLower_v<MT4> )
951  ?( IsStrictlyLower_v<MT4> ? B.lowerBound(i+2UL,j) : B.upperBound(i+2UL,j) )
952  :( B.end(j) ) );
953 
954  for( ; element!=end; ++element ) {
955  C(i ,j) += A(i ,element->index()) * element->value();
956  C(i+1UL,j) += A(i+1UL,element->index()) * element->value();
957  }
958  }
959  }
960 
961  for( ; i<M; ++i ) {
962  for( size_t j=( UPP ? i : 0UL ); j<( LOW ? i+1UL : N ); ++j )
963  {
964  auto element( ( IsUpper_v<MT4> )
965  ?( IsStrictlyUpper_v<MT4> ? B.upperBound(i,j) : B.lowerBound(i,j) )
966  :( B.begin(j) ) );
967  const auto end( ( IsLower_v<MT4> )
968  ?( IsStrictlyLower_v<MT4> ? B.lowerBound(i,j) : B.upperBound(i,j) )
969  :( B.end(j) ) );
970 
971  for( ; element!=end; ++element )
972  C(i,j) += A(i,element->index()) * element->value();
973  }
974  }
975  }
976  }
978  //**********************************************************************************************
979 
980  //**Optimized addition assignment to dense matrices*********************************************
994  template< typename MT3 // Type of the left-hand side target matrix
995  , typename MT4 // Type of the left-hand side matrix operand
996  , typename MT5 > // Type of the right-hand side matrix operand
997  static inline auto selectAddAssignKernel( MT3& C, const MT4& A, const MT5& B )
998  -> EnableIf_t< UseOptimizedKernel_v<MT3,MT4,MT5> >
999  {
1000  const size_t M( A.rows() );
1001  const size_t N( B.columns() );
1002 
1003  BLAZE_INTERNAL_ASSERT( !( LOW || UPP ) || M == N, "Broken invariant detected" );
1004 
1005  {
1006  size_t i( 0UL );
1007 
1008  for( ; (i+4UL) <= M; i+=4UL ) {
1009  for( size_t j=( UPP ? i : 0UL ); j<( LOW ? i+4UL : N ); ++j )
1010  {
1011  auto element( ( IsUpper_v<MT4> )
1012  ?( IsStrictlyUpper_v<MT4> ? B.upperBound(i,j) : B.lowerBound(i,j) )
1013  :( B.begin(j) ) );
1014  const auto end( ( IsLower_v<MT4> )
1015  ?( IsStrictlyLower_v<MT4> ? B.lowerBound(i+4UL,j) : B.upperBound(i+4UL,j) )
1016  :( B.end(j) ) );
1017 
1018  const size_t nonzeros( end - element );
1019  const size_t kpos( nonzeros & size_t(-4) );
1020  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
1021 
1022  for( size_t k=0UL; k<kpos; k+=4UL )
1023  {
1024  const size_t j1( element->index() );
1025  const ET2 v1( element->value() );
1026  ++element;
1027  const size_t j2( element->index() );
1028  const ET2 v2( element->value() );
1029  ++element;
1030  const size_t j3( element->index() );
1031  const ET2 v3( element->value() );
1032  ++element;
1033  const size_t j4( element->index() );
1034  const ET2 v4( element->value() );
1035  ++element;
1036 
1037  BLAZE_INTERNAL_ASSERT( j1 < j2 && j2 < j3 && j3 < j4, "Invalid sparse matrix index detected" );
1038 
1039  C(i ,j) += A(i ,j1) * v1 + A(i ,j2) * v2 + A(i ,j3) * v3 + A(i ,j4) * v4;
1040  C(i+1UL,j) += A(i+1UL,j1) * v1 + A(i+1UL,j2) * v2 + A(i+1UL,j3) * v3 + A(i+1UL,j4) * v4;
1041  C(i+2UL,j) += A(i+2UL,j1) * v1 + A(i+2UL,j2) * v2 + A(i+2UL,j3) * v3 + A(i+2UL,j4) * v4;
1042  C(i+3UL,j) += A(i+3UL,j1) * v1 + A(i+3UL,j2) * v2 + A(i+3UL,j3) * v3 + A(i+3UL,j4) * v4;
1043  }
1044 
1045  for( ; element!=end; ++element )
1046  {
1047  const size_t j1( element->index() );
1048  const ET2 v1( element->value() );
1049 
1050  C(i ,j) += A(i ,j1) * v1;
1051  C(i+1UL,j) += A(i+1UL,j1) * v1;
1052  C(i+2UL,j) += A(i+2UL,j1) * v1;
1053  C(i+3UL,j) += A(i+3UL,j1) * v1;
1054  }
1055  }
1056  }
1057 
1058  for( ; (i+2UL) <= M; i+=2UL ) {
1059  for( size_t j=( UPP ? i : 0UL ); j<( LOW ? i+2UL : N ); ++j )
1060  {
1061  auto element( ( IsUpper_v<MT4> )
1062  ?( IsStrictlyUpper_v<MT4> ? B.upperBound(i,j) : B.lowerBound(i,j) )
1063  :( B.begin(j) ) );
1064  const auto end( ( IsLower_v<MT4> )
1065  ?( IsStrictlyLower_v<MT4> ? B.lowerBound(i+2UL,j) : B.upperBound(i+2UL,j) )
1066  :( B.end(j) ) );
1067 
1068  const size_t nonzeros( end - element );
1069  const size_t kpos( nonzeros & size_t(-4) );
1070  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
1071 
1072  for( size_t k=0UL; k<kpos; k+=4UL )
1073  {
1074  const size_t j1( element->index() );
1075  const ET2 v1( element->value() );
1076  ++element;
1077  const size_t j2( element->index() );
1078  const ET2 v2( element->value() );
1079  ++element;
1080  const size_t j3( element->index() );
1081  const ET2 v3( element->value() );
1082  ++element;
1083  const size_t j4( element->index() );
1084  const ET2 v4( element->value() );
1085  ++element;
1086 
1087  BLAZE_INTERNAL_ASSERT( j1 < j2 && j2 < j3 && j3 < j4, "Invalid sparse matrix index detected" );
1088 
1089  C(i ,j) += A(i ,j1) * v1 + A(i ,j2) * v2 + A(i ,j3) * v3 + A(i ,j4) * v4;
1090  C(i+1UL,j) += A(i+1UL,j1) * v1 + A(i+1UL,j2) * v2 + A(i+1UL,j3) * v3 + A(i+1UL,j4) * v4;
1091  }
1092 
1093  for( ; element!=end; ++element )
1094  {
1095  const size_t j1( element->index() );
1096  const ET2 v1( element->value() );
1097 
1098  C(i ,j) += A(i ,j1) * v1;
1099  C(i+1UL,j) += A(i+1UL,j1) * v1;
1100  }
1101  }
1102  }
1103 
1104  for( ; i<M; ++i ) {
1105  for( size_t j=( UPP ? i : 0UL ); j<( LOW ? i+1UL : N ); ++j )
1106  {
1107  auto element( ( IsUpper_v<MT4> )
1108  ?( IsStrictlyUpper_v<MT4> ? B.upperBound(i,j) : B.lowerBound(i,j) )
1109  :( B.begin(j) ) );
1110  const auto end( ( IsLower_v<MT4> )
1111  ?( IsStrictlyLower_v<MT4> ? B.lowerBound(i,j) : B.upperBound(i,j) )
1112  :( B.end(j) ) );
1113 
1114  const size_t nonzeros( end - element );
1115  const size_t kpos( nonzeros & size_t(-4) );
1116  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
1117 
1118  for( size_t k=0UL; k<kpos; k+=4UL )
1119  {
1120  const size_t j1( element->index() );
1121  const ET2 v1( element->value() );
1122  ++element;
1123  const size_t j2( element->index() );
1124  const ET2 v2( element->value() );
1125  ++element;
1126  const size_t j3( element->index() );
1127  const ET2 v3( element->value() );
1128  ++element;
1129  const size_t j4( element->index() );
1130  const ET2 v4( element->value() );
1131  ++element;
1132 
1133  BLAZE_INTERNAL_ASSERT( j1 < j2 && j2 < j3 && j3 < j4, "Invalid sparse matrix index detected" );
1134 
1135  C(i,j) += A(i,j1) * v1 + A(i,j2) * v2 + A(i,j3) * v3 + A(i,j4) * v4;
1136  }
1137 
1138  for( ; element!=end; ++element )
1139  {
1140  const size_t j1( element->index() );
1141  const ET2 v1( element->value() );
1142 
1143  C(i,j) += A(i,j1) * v1;
1144  }
1145  }
1146  }
1147  }
1148  }
1150  //**********************************************************************************************
1151 
1152  //**Restructuring addition assignment***********************************************************
1167  template< typename MT // Type of the target matrix
1168  , bool SO > // Storage order of the target matrix
1169  friend inline auto addAssign( Matrix<MT,SO>& lhs, const DMatTSMatMultExpr& rhs )
1170  -> EnableIf_t< CanExploitSymmetry_v<MT,MT1,MT2> >
1171  {
1173 
1175 
1176  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1177  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1178 
1179  const ForwardFunctor fwd;
1180 
1181  addAssign( ~lhs, fwd( trans( rhs.lhs_ ) * rhs.rhs_ ) );
1182  }
1184  //**********************************************************************************************
1185 
1186  //**Addition assignment to sparse matrices******************************************************
1187  // No special implementation for the addition assignment to sparse matrices.
1188  //**********************************************************************************************
1189 
1190  //**Subtraction assignment to dense matrices****************************************************
1203  template< typename MT // Type of the target dense matrix
1204  , bool SO > // Storage order of the target dense matrix
1205  friend inline auto subAssign( DenseMatrix<MT,SO>& lhs, const DMatTSMatMultExpr& rhs )
1206  -> DisableIf_t< CanExploitSymmetry_v<MT,MT1,MT2> >
1207  {
1209 
1210  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1211  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1212 
1213  LT A( serial( rhs.lhs_ ) ); // Evaluation of the left-hand side dense matrix operand
1214  RT B( serial( rhs.rhs_ ) ); // Evaluation of the right-hand side sparse matrix operand
1215 
1216  BLAZE_INTERNAL_ASSERT( A.rows() == (~lhs).rows() , "Invalid number of rows" );
1217  BLAZE_INTERNAL_ASSERT( A.columns() == B.rows() , "Invalid matrix sizes" );
1218  BLAZE_INTERNAL_ASSERT( B.columns() == (~lhs).columns(), "Invalid number of columns" );
1219 
1220  DMatTSMatMultExpr::selectSubAssignKernel( ~lhs, A, B );
1221  }
1223  //**********************************************************************************************
1224 
1225  //**Default subtraction assignment to dense matrices********************************************
1239  template< typename MT3 // Type of the left-hand side target matrix
1240  , typename MT4 // Type of the left-hand side matrix operand
1241  , typename MT5 > // Type of the right-hand side matrix operand
1242  static inline auto selectSubAssignKernel( MT3& C, const MT4& A, const MT5& B )
1243  -> DisableIf_t< UseOptimizedKernel_v<MT3,MT4,MT5> >
1244  {
1245  const size_t M( A.rows() );
1246  const size_t N( B.columns() );
1247 
1248  BLAZE_INTERNAL_ASSERT( !( LOW || UPP ) || M == N, "Broken invariant detected" );
1249 
1250  {
1251  size_t i( 0UL );
1252 
1253  for( ; (i+4UL) <= M; i+=4UL ) {
1254  for( size_t j=( UPP ? i : 0UL ); j<( LOW ? i+4UL : N ); ++j )
1255  {
1256  auto element( ( IsUpper_v<MT4> )
1257  ?( IsStrictlyUpper_v<MT4> ? B.upperBound(i,j) : B.lowerBound(i,j) )
1258  :( B.begin(j) ) );
1259  const auto end( ( IsLower_v<MT4> )
1260  ?( IsStrictlyLower_v<MT4> ? B.lowerBound(i+4UL,j) : B.upperBound(i+4UL,j) )
1261  :( B.end(j) ) );
1262 
1263  for( ; element!=end; ++element ) {
1264  C(i ,j) -= A(i ,element->index()) * element->value();
1265  C(i+1UL,j) -= A(i+1UL,element->index()) * element->value();
1266  C(i+2UL,j) -= A(i+2UL,element->index()) * element->value();
1267  C(i+3UL,j) -= A(i+3UL,element->index()) * element->value();
1268  }
1269  }
1270  }
1271 
1272  for( ; (i+2UL) <= M; i+=2UL ) {
1273  for( size_t j=( UPP ? i : 0UL ); j<( LOW ? i+2UL : N ); ++j )
1274  {
1275  auto element( ( IsUpper_v<MT4> )
1276  ?( IsStrictlyUpper_v<MT4> ? B.upperBound(i,j) : B.lowerBound(i,j) )
1277  :( B.begin(j) ) );
1278  const auto end( ( IsLower_v<MT4> )
1279  ?( IsStrictlyLower_v<MT4> ? B.lowerBound(i+2UL,j) : B.upperBound(i+2UL,j) )
1280  :( B.end(j) ) );
1281 
1282  for( ; element!=end; ++element ) {
1283  C(i ,j) -= A(i ,element->index()) * element->value();
1284  C(i+1UL,j) -= A(i+1UL,element->index()) * element->value();
1285  }
1286  }
1287  }
1288 
1289  for( ; i<M; ++i ) {
1290  for( size_t j=( UPP ? i : 0UL ); j<( LOW ? i+1UL : N ); ++j )
1291  {
1292  auto element( ( IsUpper_v<MT4> )
1293  ?( IsStrictlyUpper_v<MT4> ? B.upperBound(i,j) : B.lowerBound(i,j) )
1294  :( B.begin(j) ) );
1295  const auto end( ( IsLower_v<MT4> )
1296  ?( IsStrictlyLower_v<MT4> ? B.lowerBound(i,j) : B.upperBound(i,j) )
1297  :( B.end(j) ) );
1298 
1299  for( ; element!=end; ++element )
1300  C(i,j) -= A(i,element->index()) * element->value();
1301  }
1302  }
1303  }
1304  }
1306  //**********************************************************************************************
1307 
1308  //**Optimized subtraction assignment to dense matrices******************************************
1322  template< typename MT3 // Type of the left-hand side target matrix
1323  , typename MT4 // Type of the left-hand side matrix operand
1324  , typename MT5 > // Type of the right-hand side matrix operand
1325  static inline auto selectSubAssignKernel( MT3& C, const MT4& A, const MT5& B )
1326  -> EnableIf_t< UseOptimizedKernel_v<MT3,MT4,MT5> >
1327  {
1328  const size_t M( A.rows() );
1329  const size_t N( B.columns() );
1330 
1331  BLAZE_INTERNAL_ASSERT( !( LOW || UPP ) || M == N, "Broken invariant detected" );
1332 
1333  {
1334  size_t i( 0UL );
1335 
1336  for( ; (i+4UL) <= M; i+=4UL ) {
1337  for( size_t j=( UPP ? i : 0UL ); j<( LOW ? i+4UL : N ); ++j )
1338  {
1339  auto element( ( IsUpper_v<MT4> )
1340  ?( IsStrictlyUpper_v<MT4> ? B.upperBound(i,j) : B.lowerBound(i,j) )
1341  :( B.begin(j) ) );
1342  const auto end( ( IsLower_v<MT4> )
1343  ?( IsStrictlyLower_v<MT4> ? B.lowerBound(i+4UL,j) : B.upperBound(i+4UL,j) )
1344  :( B.end(j) ) );
1345 
1346  const size_t nonzeros( end - element );
1347  const size_t kpos( nonzeros & size_t(-4) );
1348  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
1349 
1350  for( size_t k=0UL; k<kpos; k+=4UL )
1351  {
1352  const size_t j1( element->index() );
1353  const ET2 v1( element->value() );
1354  ++element;
1355  const size_t j2( element->index() );
1356  const ET2 v2( element->value() );
1357  ++element;
1358  const size_t j3( element->index() );
1359  const ET2 v3( element->value() );
1360  ++element;
1361  const size_t j4( element->index() );
1362  const ET2 v4( element->value() );
1363  ++element;
1364 
1365  BLAZE_INTERNAL_ASSERT( j1 < j2 && j2 < j3 && j3 < j4, "Invalid sparse matrix index detected" );
1366 
1367  C(i ,j) -= A(i ,j1) * v1 + A(i ,j2) * v2 + A(i ,j3) * v3 + A(i ,j4) * v4;
1368  C(i+1UL,j) -= A(i+1UL,j1) * v1 + A(i+1UL,j2) * v2 + A(i+1UL,j3) * v3 + A(i+1UL,j4) * v4;
1369  C(i+2UL,j) -= A(i+2UL,j1) * v1 + A(i+2UL,j2) * v2 + A(i+2UL,j3) * v3 + A(i+2UL,j4) * v4;
1370  C(i+3UL,j) -= A(i+3UL,j1) * v1 + A(i+3UL,j2) * v2 + A(i+3UL,j3) * v3 + A(i+3UL,j4) * v4;
1371  }
1372 
1373  for( ; element!=end; ++element )
1374  {
1375  const size_t j1( element->index() );
1376  const ET2 v1( element->value() );
1377 
1378  C(i ,j) -= A(i ,j1) * v1;
1379  C(i+1UL,j) -= A(i+1UL,j1) * v1;
1380  C(i+2UL,j) -= A(i+2UL,j1) * v1;
1381  C(i+3UL,j) -= A(i+3UL,j1) * v1;
1382  }
1383  }
1384  }
1385 
1386  for( ; (i+2UL) <= M; i+=2UL ) {
1387  for( size_t j=( UPP ? i : 0UL ); j<( LOW ? i+2UL : N ); ++j )
1388  {
1389  auto element( ( IsUpper_v<MT4> )
1390  ?( IsStrictlyUpper_v<MT4> ? B.upperBound(i,j) : B.lowerBound(i,j) )
1391  :( B.begin(j) ) );
1392  const auto end( ( IsLower_v<MT4> )
1393  ?( IsStrictlyLower_v<MT4> ? B.lowerBound(i+2UL,j) : B.upperBound(i+2UL,j) )
1394  :( B.end(j) ) );
1395 
1396  const size_t nonzeros( end - element );
1397  const size_t kpos( nonzeros & size_t(-4) );
1398  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
1399 
1400  for( size_t k=0UL; k<kpos; k+=4UL )
1401  {
1402  const size_t j1( element->index() );
1403  const ET2 v1( element->value() );
1404  ++element;
1405  const size_t j2( element->index() );
1406  const ET2 v2( element->value() );
1407  ++element;
1408  const size_t j3( element->index() );
1409  const ET2 v3( element->value() );
1410  ++element;
1411  const size_t j4( element->index() );
1412  const ET2 v4( element->value() );
1413  ++element;
1414 
1415  BLAZE_INTERNAL_ASSERT( j1 < j2 && j2 < j3 && j3 < j4, "Invalid sparse matrix index detected" );
1416 
1417  C(i ,j) -= A(i ,j1) * v1 + A(i ,j2) * v2 + A(i ,j3) * v3 + A(i ,j4) * v4;
1418  C(i+1UL,j) -= A(i+1UL,j1) * v1 + A(i+1UL,j2) * v2 + A(i+1UL,j3) * v3 + A(i+1UL,j4) * v4;
1419  }
1420 
1421  for( ; element!=end; ++element )
1422  {
1423  const size_t j1( element->index() );
1424  const ET2 v1( element->value() );
1425 
1426  C(i ,j) -= A(i ,j1) * v1;
1427  C(i+1UL,j) -= A(i+1UL,j1) * v1;
1428  }
1429  }
1430  }
1431 
1432  for( ; i<M; ++i ) {
1433  for( size_t j=( UPP ? i : 0UL ); j<( LOW ? i+1UL : N ); ++j )
1434  {
1435  auto element( ( IsUpper_v<MT4> )
1436  ?( IsStrictlyUpper_v<MT4> ? B.upperBound(i,j) : B.lowerBound(i,j) )
1437  :( B.begin(j) ) );
1438  const auto end( ( IsLower_v<MT4> )
1439  ?( IsStrictlyLower_v<MT4> ? B.lowerBound(i,j) : B.upperBound(i,j) )
1440  :( B.end(j) ) );
1441 
1442  const size_t nonzeros( end - element );
1443  const size_t kpos( nonzeros & size_t(-4) );
1444  BLAZE_INTERNAL_ASSERT( ( nonzeros - ( nonzeros % 4UL ) ) == kpos, "Invalid end calculation" );
1445 
1446  for( size_t k=0UL; k<kpos; k+=4UL )
1447  {
1448  const size_t j1( element->index() );
1449  const ET2 v1( element->value() );
1450  ++element;
1451  const size_t j2( element->index() );
1452  const ET2 v2( element->value() );
1453  ++element;
1454  const size_t j3( element->index() );
1455  const ET2 v3( element->value() );
1456  ++element;
1457  const size_t j4( element->index() );
1458  const ET2 v4( element->value() );
1459  ++element;
1460 
1461  BLAZE_INTERNAL_ASSERT( j1 < j2 && j2 < j3 && j3 < j4, "Invalid sparse matrix index detected" );
1462 
1463  C(i,j) -= A(i,j1) * v1 + A(i,j2) * v2 + A(i,j3) * v3 + A(i,j4) * v4;
1464  }
1465 
1466  for( ; element!=end; ++element )
1467  {
1468  const size_t j1( element->index() );
1469  const ET2 v1( element->value() );
1470 
1471  C(i,j) -= A(i,j1) * v1;
1472  }
1473  }
1474  }
1475  }
1476  }
1478  //**********************************************************************************************
1479 
1480  //**Restructuring subtraction assignment********************************************************
1495  template< typename MT // Type of the target matrix
1496  , bool SO > // Storage order of the target matrix
1497  friend inline auto subAssign( Matrix<MT,SO>& lhs, const DMatTSMatMultExpr& rhs )
1498  -> EnableIf_t< CanExploitSymmetry_v<MT,MT1,MT2> >
1499  {
1501 
1503 
1504  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1505  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1506 
1507  const ForwardFunctor fwd;
1508 
1509  subAssign( ~lhs, fwd( trans( rhs.lhs_ ) * rhs.rhs_ ) );
1510  }
1512  //**********************************************************************************************
1513 
1514  //**Subtraction assignment to sparse matrices***************************************************
1515  // No special implementation for the subtraction assignment to sparse matrices.
1516  //**********************************************************************************************
1517 
1518  //**Schur product assignment to dense matrices**************************************************
1531  template< typename MT // Type of the target dense matrix
1532  , bool SO > // Storage order of the target dense matrix
1533  friend inline void schurAssign( DenseMatrix<MT,SO>& lhs, const DMatTSMatMultExpr& rhs )
1534  {
1536 
1540 
1541  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1542  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1543 
1544  const ResultType tmp( serial( rhs ) );
1545  schurAssign( ~lhs, tmp );
1546  }
1548  //**********************************************************************************************
1549 
1550  //**Schur product assignment to sparse matrices*************************************************
1551  // No special implementation for the Schur product assignment to sparse matrices.
1552  //**********************************************************************************************
1553 
1554  //**Multiplication assignment to dense matrices*************************************************
1555  // No special implementation for the multiplication assignment to dense matrices.
1556  //**********************************************************************************************
1557 
1558  //**Multiplication assignment to sparse matrices************************************************
1559  // No special implementation for the multiplication assignment to sparse matrices.
1560  //**********************************************************************************************
1561 
1562  //**SMP assignment to dense matrices************************************************************
1577  template< typename MT // Type of the target dense matrix
1578  , bool SO > // Storage order of the target dense matrix
1579  friend inline auto smpAssign( DenseMatrix<MT,SO>& lhs, const DMatTSMatMultExpr& rhs )
1580  -> EnableIf_t< IsEvaluationRequired_v<MT,MT1,MT2> >
1581  {
1583 
1584  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1585  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1586 
1587  LT A( rhs.lhs_ ); // Evaluation of the left-hand side dense matrix operand
1588  RT B( rhs.rhs_ ); // Evaluation of the right-hand side sparse matrix operand
1589 
1590  BLAZE_INTERNAL_ASSERT( A.rows() == (~lhs).rows() , "Invalid number of rows" );
1591  BLAZE_INTERNAL_ASSERT( A.columns() == B.rows() , "Invalid matrix sizes" );
1592  BLAZE_INTERNAL_ASSERT( B.columns() == (~lhs).columns(), "Invalid number of columns" );
1593 
1594  smpAssign( ~lhs, A * B );
1595  }
1597  //**********************************************************************************************
1598 
1599  //**SMP assignment to sparse matrices***********************************************************
1614  template< typename MT // Type of the target sparse matrix
1615  , bool SO > // Storage order of the target sparse matrix
1616  friend inline auto smpAssign( SparseMatrix<MT,SO>& lhs, const DMatTSMatMultExpr& rhs )
1617  -> EnableIf_t< IsEvaluationRequired_v<MT,MT1,MT2> >
1618  {
1620 
1621  using TmpType = If_t< SO, OppositeType, ResultType >;
1622 
1629 
1630  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1631  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1632 
1633  const ForwardFunctor fwd;
1634 
1635  const TmpType tmp( rhs );
1636  smpAssign( ~lhs, fwd( tmp ) );
1637  }
1639  //**********************************************************************************************
1640 
1641  //**Restructuring SMP assignment****************************************************************
1656  template< typename MT // Type of the target matrix
1657  , bool SO > // Storage order of the target matrix
1658  friend inline auto smpAssign( Matrix<MT,SO>& lhs, const DMatTSMatMultExpr& rhs )
1659  -> EnableIf_t< CanExploitSymmetry_v<MT,MT1,MT2> >
1660  {
1662 
1664 
1665  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1666  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1667 
1668  const ForwardFunctor fwd;
1669 
1670  smpAssign( ~lhs, fwd( trans( rhs.lhs_ ) * rhs.rhs_ ) );
1671  }
1673  //**********************************************************************************************
1674 
1675  //**SMP addition assignment to dense matrices***************************************************
1691  template< typename MT // Type of the target dense matrix
1692  , bool SO > // Storage order of the target dense matrix
1693  friend inline auto smpAddAssign( DenseMatrix<MT,SO>& lhs, const DMatTSMatMultExpr& rhs )
1694  -> EnableIf_t< IsEvaluationRequired_v<MT,MT1,MT2> >
1695  {
1697 
1698  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1699  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1700 
1701  LT A( rhs.lhs_ ); // Evaluation of the left-hand side dense matrix operand
1702  RT B( rhs.rhs_ ); // Evaluation of the right-hand side sparse matrix operand
1703 
1704  BLAZE_INTERNAL_ASSERT( A.rows() == (~lhs).rows() , "Invalid number of rows" );
1705  BLAZE_INTERNAL_ASSERT( A.columns() == B.rows() , "Invalid matrix sizes" );
1706  BLAZE_INTERNAL_ASSERT( B.columns() == (~lhs).columns(), "Invalid number of columns" );
1707 
1708  smpAddAssign( ~lhs, A * B );
1709  }
1711  //**********************************************************************************************
1712 
1713  //**Restructuring SMP addition assignment*******************************************************
1728  template< typename MT // Type of the target matrix
1729  , bool SO > // Storage order of the target matrix
1730  friend inline auto smpAddAssign( Matrix<MT,SO>& lhs, const DMatTSMatMultExpr& rhs )
1731  -> EnableIf_t< CanExploitSymmetry_v<MT,MT1,MT2> >
1732  {
1734 
1736 
1737  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1738  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1739 
1740  const ForwardFunctor fwd;
1741 
1742  smpAddAssign( ~lhs, fwd( trans( rhs.lhs_ ) * rhs.rhs_ ) );
1743  }
1745  //**********************************************************************************************
1746 
1747  //**SMP addition assignment to sparse matrices**************************************************
1748  // No special implementation for the SMP addition assignment to sparse matrices.
1749  //**********************************************************************************************
1750 
1751  //**SMP subtraction assignment to dense matrices************************************************
1767  template< typename MT // Type of the target dense matrix
1768  , bool SO > // Storage order of the target dense matrix
1769  friend inline auto smpSubAssign( DenseMatrix<MT,SO>& lhs, const DMatTSMatMultExpr& rhs )
1770  -> EnableIf_t< IsEvaluationRequired_v<MT,MT1,MT2> >
1771  {
1773 
1774  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1775  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1776 
1777  LT A( rhs.lhs_ ); // Evaluation of the left-hand side dense matrix operand
1778  RT B( rhs.rhs_ ); // Evaluation of the right-hand side sparse matrix operand
1779 
1780  BLAZE_INTERNAL_ASSERT( A.rows() == (~lhs).rows() , "Invalid number of rows" );
1781  BLAZE_INTERNAL_ASSERT( A.columns() == B.rows() , "Invalid matrix sizes" );
1782  BLAZE_INTERNAL_ASSERT( B.columns() == (~lhs).columns(), "Invalid number of columns" );
1783 
1784  smpSubAssign( ~lhs, A * B );
1785  }
1787  //**********************************************************************************************
1788 
1789  //**Restructuring SMP subtraction assignment****************************************************
1804  template< typename MT // Type of the target matrix
1805  , bool SO > // Storage order of the target matrix
1806  friend inline auto smpSubAssign( Matrix<MT,SO>& lhs, const DMatTSMatMultExpr& rhs )
1807  -> EnableIf_t< CanExploitSymmetry_v<MT,MT1,MT2> >
1808  {
1810 
1812 
1813  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1814  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1815 
1816  const ForwardFunctor fwd;
1817 
1818  smpSubAssign( ~lhs, fwd( trans( rhs.lhs_ ) * rhs.rhs_ ) );
1819  }
1821  //**********************************************************************************************
1822 
1823  //**SMP subtraction assignment to sparse matrices***********************************************
1824  // No special implementation for the SMP subtraction assignment to sparse matrices.
1825  //**********************************************************************************************
1826 
1827  //**SMP Schur product assignment to dense matrices**********************************************
1840  template< typename MT // Type of the target dense matrix
1841  , bool SO > // Storage order of the target dense matrix
1842  friend inline void smpSchurAssign( DenseMatrix<MT,SO>& lhs, const DMatTSMatMultExpr& rhs )
1843  {
1845 
1849 
1850  BLAZE_INTERNAL_ASSERT( (~lhs).rows() == rhs.rows() , "Invalid number of rows" );
1851  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == rhs.columns(), "Invalid number of columns" );
1852 
1853  const ResultType tmp( rhs );
1854  smpSchurAssign( ~lhs, tmp );
1855  }
1857  //**********************************************************************************************
1858 
1859  //**SMP Schur product assignment to sparse matrices*********************************************
1860  // No special implementation for the SMP Schur product assignment to sparse matrices.
1861  //**********************************************************************************************
1862 
1863  //**SMP multiplication assignment to dense matrices*********************************************
1864  // No special implementation for the SMP multiplication assignment to dense matrices.
1865  //**********************************************************************************************
1866 
1867  //**SMP multiplication assignment to sparse matrices********************************************
1868  // No special implementation for the SMP multiplication assignment to sparse matrices.
1869  //**********************************************************************************************
1870 
1871  //**Compile time checks*************************************************************************
1880  //**********************************************************************************************
1881 };
1882 //*************************************************************************************************
1883 
1884 
1885 
1886 
1887 //=================================================================================================
1888 //
1889 // GLOBAL BINARY ARITHMETIC OPERATORS
1890 //
1891 //=================================================================================================
1892 
1893 //*************************************************************************************************
1906 template< typename MT1 // Type of the left-hand side dense matrix
1907  , typename MT2 // Type of the right-hand side sparse matrix
1908  , DisableIf_t< ( IsIdentity_v<MT2> &&
1909  IsSame_v< ElementType_t<MT1>, ElementType_t<MT2> > ) ||
1910  IsZero_v<MT2> >* = nullptr >
1911 inline const DMatTSMatMultExpr<MT1,MT2,false,false,false,false>
1912  dmattsmatmult( const DenseMatrix<MT1,false>& lhs, const SparseMatrix<MT2,true>& rhs )
1913 {
1915 
1916  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == (~rhs).rows(), "Invalid matrix sizes" );
1917 
1918  return DMatTSMatMultExpr<MT1,MT2,false,false,false,false>( ~lhs, ~rhs );
1919 }
1921 //*************************************************************************************************
1922 
1923 
1924 //*************************************************************************************************
1938 template< typename MT1 // Type of the left-hand side dense matrix
1939  , typename MT2 // Type of the right-hand side sparse matrix
1940  , EnableIf_t< IsIdentity_v<MT2> &&
1941  IsSame_v< ElementType_t<MT1>, ElementType_t<MT2> > >* = nullptr >
1942 inline const MT1&
1943  dmattsmatmult( const DenseMatrix<MT1,false>& lhs, const SparseMatrix<MT2,true>& rhs )
1944 {
1946 
1947  UNUSED_PARAMETER( rhs );
1948 
1949  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == (~rhs).rows(), "Invalid matrix sizes" );
1950 
1951  return (~lhs);
1952 }
1954 //*************************************************************************************************
1955 
1956 
1957 //*************************************************************************************************
1970 template< typename MT1 // Type of the left-hand side dense matrix
1971  , typename MT2 // Type of the right-hand side sparse matrix
1972  , EnableIf_t< IsZero_v<MT2> >* = nullptr >
1973 inline decltype(auto)
1974  dmattsmatmult( const DenseMatrix<MT1,false>& lhs, const SparseMatrix<MT2,true>& rhs )
1975 {
1977 
1978  BLAZE_INTERNAL_ASSERT( (~lhs).columns() == (~rhs).rows(), "Invalid matrix sizes" );
1979 
1980  using ReturnType = const MultTrait_t< ResultType_t<MT1>, ResultType_t<MT2> >;
1981 
1984 
1985  return ReturnType( (~lhs).rows(), (~rhs).columns() );
1986 }
1988 //*************************************************************************************************
1989 
1990 
1991 //*************************************************************************************************
2021 template< typename MT1 // Type of the left-hand side dense matrix
2022  , typename MT2 > // Type of the right-hand side sparse matrix
2023 inline decltype(auto)
2024  operator*( const DenseMatrix<MT1,false>& lhs, const SparseMatrix<MT2,true>& rhs )
2025 {
2027 
2028  if( (~lhs).columns() != (~rhs).rows() ) {
2029  BLAZE_THROW_INVALID_ARGUMENT( "Matrix sizes do not match" );
2030  }
2031 
2032  return dmattsmatmult( ~lhs, ~rhs );
2033 }
2034 //*************************************************************************************************
2035 
2036 
2037 
2038 
2039 //=================================================================================================
2040 //
2041 // GLOBAL FUNCTIONS
2042 //
2043 //=================================================================================================
2044 
2045 //*************************************************************************************************
2069 template< typename MT1 // Type of the left-hand side dense matrix
2070  , typename MT2 // Type of the right-hand side sparse matrix
2071  , bool SF // Symmetry flag
2072  , bool HF // Hermitian flag
2073  , bool LF // Lower flag
2074  , bool UF > // Upper flag
2075 inline decltype(auto) declsym( const DMatTSMatMultExpr<MT1,MT2,SF,HF,LF,UF>& dm )
2076 {
2078 
2079  if( !isSquare( dm ) ) {
2080  BLAZE_THROW_INVALID_ARGUMENT( "Invalid symmetric matrix specification" );
2081  }
2082 
2083  using ReturnType = const DMatTSMatMultExpr<MT1,MT2,true,HF,LF,UF>;
2084  return ReturnType( dm.leftOperand(), dm.rightOperand() );
2085 }
2087 //*************************************************************************************************
2088 
2089 
2090 //*************************************************************************************************
2114 template< typename MT1 // Type of the left-hand side dense matrix
2115  , typename MT2 // Type of the right-hand side sparse matrix
2116  , bool SF // Symmetry flag
2117  , bool HF // Hermitian flag
2118  , bool LF // Lower flag
2119  , bool UF > // Upper flag
2120 inline decltype(auto) declherm( const DMatTSMatMultExpr<MT1,MT2,SF,HF,LF,UF>& dm )
2121 {
2123 
2124  if( !isSquare( dm ) ) {
2125  BLAZE_THROW_INVALID_ARGUMENT( "Invalid Hermitian matrix specification" );
2126  }
2127 
2128  using ReturnType = const DMatTSMatMultExpr<MT1,MT2,SF,true,LF,UF>;
2129  return ReturnType( dm.leftOperand(), dm.rightOperand() );
2130 }
2132 //*************************************************************************************************
2133 
2134 
2135 //*************************************************************************************************
2159 template< typename MT1 // Type of the left-hand side dense matrix
2160  , typename MT2 // Type of the right-hand side sparse matrix
2161  , bool SF // Symmetry flag
2162  , bool HF // Hermitian flag
2163  , bool LF // Lower flag
2164  , bool UF > // Upper flag
2165 inline decltype(auto) decllow( const DMatTSMatMultExpr<MT1,MT2,SF,HF,LF,UF>& dm )
2166 {
2168 
2169  if( !isSquare( dm ) ) {
2170  BLAZE_THROW_INVALID_ARGUMENT( "Invalid lower matrix specification" );
2171  }
2172 
2173  using ReturnType = const DMatTSMatMultExpr<MT1,MT2,SF,HF,true,UF>;
2174  return ReturnType( dm.leftOperand(), dm.rightOperand() );
2175 }
2177 //*************************************************************************************************
2178 
2179 
2180 //*************************************************************************************************
2204 template< typename MT1 // Type of the left-hand side dense matrix
2205  , typename MT2 // Type of the right-hand side sparse matrix
2206  , bool SF // Symmetry flag
2207  , bool HF // Hermitian flag
2208  , bool LF // Lower flag
2209  , bool UF > // Upper flag
2210 inline decltype(auto) declupp( const DMatTSMatMultExpr<MT1,MT2,SF,HF,LF,UF>& dm )
2211 {
2213 
2214  if( !isSquare( dm ) ) {
2215  BLAZE_THROW_INVALID_ARGUMENT( "Invalid upper matrix specification" );
2216  }
2217 
2218  using ReturnType = const DMatTSMatMultExpr<MT1,MT2,SF,HF,LF,true>;
2219  return ReturnType( dm.leftOperand(), dm.rightOperand() );
2220 }
2222 //*************************************************************************************************
2223 
2224 
2225 //*************************************************************************************************
2249 template< typename MT1 // Type of the left-hand side dense matrix
2250  , typename MT2 // Type of the right-hand side sparse matrix
2251  , bool SF // Symmetry flag
2252  , bool HF // Hermitian flag
2253  , bool LF // Lower flag
2254  , bool UF > // Upper flag
2255 inline decltype(auto) decldiag( const DMatTSMatMultExpr<MT1,MT2,SF,HF,LF,UF>& dm )
2256 {
2258 
2259  if( !isSquare( dm ) ) {
2260  BLAZE_THROW_INVALID_ARGUMENT( "Invalid diagonal matrix specification" );
2261  }
2262 
2263  using ReturnType = const DMatTSMatMultExpr<MT1,MT2,SF,HF,true,true>;
2264  return ReturnType( dm.leftOperand(), dm.rightOperand() );
2265 }
2267 //*************************************************************************************************
2268 
2269 
2270 
2271 
2272 //=================================================================================================
2273 //
2274 // SIZE SPECIALIZATIONS
2275 //
2276 //=================================================================================================
2277 
2278 //*************************************************************************************************
2280 template< typename MT1, typename MT2, bool SF, bool HF, bool LF, bool UF >
2281 struct Size< DMatTSMatMultExpr<MT1,MT2,SF,HF,LF,UF>, 0UL >
2282  : public Size<MT1,0UL>
2283 {};
2284 
2285 template< typename MT1, typename MT2, bool SF, bool HF, bool LF, bool UF >
2286 struct Size< DMatTSMatMultExpr<MT1,MT2,SF,HF,LF,UF>, 1UL >
2287  : public Size<MT2,1UL>
2288 {};
2290 //*************************************************************************************************
2291 
2292 
2293 
2294 
2295 //=================================================================================================
2296 //
2297 // ISALIGNED SPECIALIZATIONS
2298 //
2299 //=================================================================================================
2300 
2301 //*************************************************************************************************
2303 template< typename MT1, typename MT2, bool SF, bool HF, bool LF, bool UF >
2304 struct IsAligned< DMatTSMatMultExpr<MT1,MT2,SF,HF,LF,UF> >
2305  : public IsAligned<MT1>
2306 {};
2308 //*************************************************************************************************
2309 
2310 } // namespace blaze
2311 
2312 #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 the UNUSED_PARAMETER function template.
Header file for basic type definitions.
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
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
size_t columns() const noexcept
Returns the current number of columns of the matrix.
Definition: DMatTSMatMultExpr.h:360
Header file for the serial shim.
Header file for the IsDiagonal type trait.
Base template for the DeclUppTrait class.
Definition: DeclUppTrait.h:134
If_t< evaluateRight, const RT2, CT2 > RT
Type for the assignment of the right-hand side sparse matrix operand.
Definition: DMatTSMatMultExpr.h:252
#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
ElementType_t< RT1 > ET1
Element type of the left-hand side dense matrix expression.
Definition: DMatTSMatMultExpr.h:135
Header file for the DeclUpp functor.
Header file for the IsSame and IsStrictlySame type traits.
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
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
DMatTSMatMultExpr(const MT1 &lhs, const MT2 &rhs) noexcept
Constructor for the DMatTSMatMultExpr class.
Definition: DMatTSMatMultExpr.h:270
LeftOperand leftOperand() const noexcept
Returns the left-hand side dense matrix operand.
Definition: DMatTSMatMultExpr.h:370
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.
ReturnType at(size_t i, size_t j) const
Checked access to the matrix elements.
Definition: DMatTSMatMultExpr.h:334
static constexpr bool evaluateLeft
Compilation switch for the composite type of the left-hand side dense matrix expression.
Definition: DMatTSMatMultExpr.h:143
Header file for the IsIdentity type trait.
bool canSMPAssign() const noexcept
Returns whether the expression can be used in SMP assignments.
Definition: DMatTSMatMultExpr.h:424
decltype(auto) declupp(const DenseMatrix< MT, SO > &dm)
Declares the given dense matrix expression dm as upper.
Definition: DMatDeclUppExpr.h:1002
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: DMatTSMatMultExpr.h:234
ResultType_t< MT1 > RT1
Result type of the left-hand side dense matrix expression.
Definition: DMatTSMatMultExpr.h:133
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.
ElementType_t< ResultType > ElementType
Resulting element type.
Definition: DMatTSMatMultExpr.h:238
Header file for the RequiresEvaluation type trait.
System settings for performance optimizations.
constexpr void UNUSED_PARAMETER(const Args &...)
Suppression of unused parameter warnings.
Definition: Unused.h:81
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
TransposeType_t< ResultType > TransposeType
Transpose type for expression template evaluations.
Definition: DMatTSMatMultExpr.h:237
ResultType_t< MT2 > RT2
Result type of the right-hand side sparse matrix expression.
Definition: DMatTSMatMultExpr.h:134
OppositeType_t< ResultType > OppositeType
Result type with opposite storage order for expression template evaluations.
Definition: DMatTSMatMultExpr.h:236
static constexpr bool simdEnabled
Compilation switch for the expression template evaluation strategy.
Definition: DMatTSMatMultExpr.h:257
Constraint on the data type.
Constraint on the data type.
RightOperand rhs_
Right-hand side sparse matrix of the multiplication expression.
Definition: DMatTSMatMultExpr.h:432
If_t< evaluateLeft, const RT1, CT1 > LT
Type for the assignment of the left-hand side dense matrix operand.
Definition: DMatTSMatMultExpr.h:249
Headerfile for the generic max algorithm.
Header file for the DisableIf class template.
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.
bool isAliased(const T *alias) const noexcept
Returns whether the expression is aliased with the given address alias.
Definition: DMatTSMatMultExpr.h:404
Header file for the If class template.
bool canAlias(const T *alias) const noexcept
Returns whether the expression can alias with the given address alias.
Definition: DMatTSMatMultExpr.h:392
#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
#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
Generic wrapper for the decllow() function.
Definition: DeclLow.h:58
static constexpr bool HERM
Flag for Hermitian matrices.
Definition: DMatTSMatMultExpr.h:153
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
static constexpr bool smpAssignable
Compilation switch for the expression template assignment strategy.
Definition: DMatTSMatMultExpr.h:260
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.
RightOperand rightOperand() const noexcept
Returns the right-hand side transpose sparse matrix operand.
Definition: DMatTSMatMultExpr.h:380
decltype(auto) decllow(const DenseMatrix< MT, SO > &dm)
Declares the given dense matrix expression dm as lower.
Definition: DMatDeclLowExpr.h:1002
CompositeType_t< MT1 > CT1
Composite type of the left-hand side dense matrix expression.
Definition: DMatTSMatMultExpr.h:137
Header file for the IsLower type trait.
Header file for the IsAligned type trait.
static constexpr bool SYM
Flag for symmetric matrices.
Definition: DMatTSMatMultExpr.h:152
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
size_t rows() const noexcept
Returns the current number of rows of the matrix.
Definition: DMatTSMatMultExpr.h:350
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.
ElementType_t< RT2 > ET2
Element type of the right-hand side sparse matrix expression.
Definition: DMatTSMatMultExpr.h:136
Constraint on the data type.
Header file for all forward declarations for expression class templates.
bool isAligned() const noexcept
Returns whether the operands of the expression are properly aligned in memory.
Definition: DMatTSMatMultExpr.h:414
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.
If_t< IsExpression_v< MT2 >, const MT2, const MT2 &> RightOperand
Composite type of the right-hand side sparse matrix expression.
Definition: DMatTSMatMultExpr.h:246
Base class for all matrix/matrix multiplication expression templates.The MatMatMultExpr class serves ...
Definition: MatMatMultExpr.h:67
static constexpr bool LOW
Flag for lower matrices.
Definition: DMatTSMatMultExpr.h:154
#define BLAZE_CONSTRAINT_MUST_NOT_BE_SYMMETRIC_MATRIX_TYPE(T)
Constraint on the data type.In case the given data type T is a symmetric matrix type, a compilation error is created.
Definition: Symmetric.h:79
#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
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
decltype(auto) row(Matrix< MT, SO > &, RRAs...)
Creating a view on a specific row of the given matrix.
Definition: Row.h:133
LeftOperand lhs_
Left-hand side dense matrix of the multiplication expression.
Definition: DMatTSMatMultExpr.h:431
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
Expression object for dense matrix-transpose sparse matrix multiplications.The DMatTSMatMultExpr clas...
Definition: DMatTSMatMultExpr.h:127
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
Constraint on the data type.
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
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
CompositeType_t< MT2 > CT2
Composite type of the right-hand side sparse matrix expression.
Definition: DMatTSMatMultExpr.h:138
#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
const ElementType ReturnType
Return type for expression template evaluations.
Definition: DMatTSMatMultExpr.h:239
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
const ResultType CompositeType
Data type for composite expression templates.
Definition: DMatTSMatMultExpr.h:240
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
static constexpr bool UPP
Flag for upper matrices.
Definition: DMatTSMatMultExpr.h:155
Header file for the IsComputation type trait class.
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
Base class for all compute expression templates.The Computation class serves as a tag for all computa...
Definition: Computation.h:66
Header file for the IntegralConstant class template.
Generic wrapper for the decldiag() function.
Definition: DeclDiag.h:58
If_t< IsExpression_v< MT1 >, const MT1, const MT1 &> LeftOperand
Composite type of the left-hand side dense matrix expression.
Definition: DMatTSMatMultExpr.h:243
Header file for the DeclHerm functor.
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.
static constexpr bool evaluateRight
Compilation switch for the composite type of the right-hand side sparse matrix expression.
Definition: DMatTSMatMultExpr.h:148
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
ReturnType operator()(size_t i, size_t j) const
2D-access to the matrix elements.
Definition: DMatTSMatMultExpr.h:285
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.