LvArray
genericTensorOps.hpp
Go to the documentation of this file.
1 /*
2  * Copyright (c) 2021, Lawrence Livermore National Security, LLC and LvArray contributors.
3  * All rights reserved.
4  * See the LICENSE file for details.
5  * SPDX-License-Identifier: (BSD-3-Clause)
6  */
7 
13 #pragma once
14 
15 // Source includes
16 #include "ArraySlice.hpp"
17 #include "math.hpp"
18 
43 
49 #define LVARRAY_TENSOROPS_INIT_LOCAL_2( EXP ) { EXP[ 0 ], EXP[ 1 ] }
50 
57 #define LVARRAY_TENSOROPS_ASSIGN_2( var, exp0, exp1 ) \
58  var[ 0 ] = exp0; var[ 1 ] = exp1
59 
64 #define LVARRAY_TENSOROPS_INIT_LOCAL_3( EXP ) { EXP[ 0 ], EXP[ 1 ], EXP[ 2 ] }
65 
73 #define LVARRAY_TENSOROPS_ASSIGN_3( var, exp0, exp1, exp2 ) \
74  var[ 0 ] = exp0; var[ 1 ] = exp1; var[ 2 ] = exp2
75 
80 #define LVARRAY_TENSOROPS_INIT_LOCAL_6( EXP ) { EXP[ 0 ], EXP[ 1 ], EXP[ 2 ], EXP[ 3 ], EXP[ 4 ], EXP[ 5 ] }
81 
92 #define LVARRAY_TENSOROPS_ASSIGN_6( var, exp0, exp1, exp2, exp3, exp4, exp5 ) \
93  var[ 0 ] = exp0; var[ 1 ] = exp1; var[ 2 ] = exp2; var[ 3 ] = exp3; var[ 4 ] = exp4; var[ 5 ] = exp5
94 
99 #define LVARRAY_TENSOROPS_INIT_LOCAL_2x2( EXP ) \
100  { { EXP[ 0 ][ 0 ], EXP[ 0 ][ 1 ] }, \
101  { EXP[ 1 ][ 0 ], EXP[ 1 ][ 1 ] } }
102 
111 #define LVARRAY_TENSOROPS_ASSIGN_2x2( var, exp00, exp01, exp10, exp11 ) \
112  var[ 0 ][ 0 ] = exp00; var[ 0 ][ 1 ] = exp01; \
113  var[ 1 ][ 0 ] = exp10; var[ 1 ][ 1 ] = exp11
114 
119 #define LVARRAY_TENSOROPS_INIT_LOCAL_3x3( EXP ) \
120  { { EXP[ 0 ][ 0 ], EXP[ 0 ][ 1 ], EXP[ 0 ][ 2 ] }, \
121  { EXP[ 1 ][ 0 ], EXP[ 1 ][ 1 ], EXP[ 1 ][ 2 ] }, \
122  { EXP[ 2 ][ 0 ], EXP[ 2 ][ 1 ], EXP[ 2 ][ 2 ] } }
123 
137 #define LVARRAY_TENSOROPS_ASSIGN_3x3( var, exp00, exp01, exp02, exp10, exp11, exp12, exp20, exp21, exp22 ) \
138  var[ 0 ][ 0 ] = exp00; var[ 0 ][ 1 ] = exp01; var[ 0 ][ 2 ] = exp02; \
139  var[ 1 ][ 0 ] = exp10; var[ 1 ][ 1 ] = exp11; var[ 1 ][ 2 ] = exp12; \
140  var[ 2 ][ 0 ] = exp20; var[ 2 ][ 1 ] = exp21; var[ 2 ][ 2 ] = exp22
141 
143 
144 
145 namespace LvArray
146 {
147 namespace tensorOps
148 {
149 namespace internal
150 {
151 
157 HAS_STATIC_MEMBER( SIZE );
158 
164 HAS_STATIC_MEMBER( NDIM );
165 
175 template< std::ptrdiff_t ISIZE, typename T >
176 LVARRAY_HOST_DEVICE inline constexpr
177 std::enable_if_t< HasStaticMember_SIZE< T > >
178 checkSizes( T const & src )
179 {
180  static_assert( ISIZE == T::SIZE,
181  "Expected the first dimension of size ISIZE, got an type of size T::SIZE." );
182  LVARRAY_UNUSED_VARIABLE( src );
183 }
184 
192 template< std::ptrdiff_t PROVIDED_SIZE, typename T, std::ptrdiff_t INFERRED_SIZE >
193 LVARRAY_HOST_DEVICE inline constexpr
194 void checkSizes( T const ( &src )[ INFERRED_SIZE ] )
195 {
196  static_assert( PROVIDED_SIZE == INFERRED_SIZE,
197  "Expected the first dimension of size PROVIDED_N, got an array of size INFERRED_N." );
198  LVARRAY_UNUSED_VARIABLE( src );
199 }
200 
210 template< std::ptrdiff_t PROVIDED_M,
211  std::ptrdiff_t PROVIDED_N,
212  typename T,
213  std::ptrdiff_t INFERRED_M,
214  std::ptrdiff_t INFERRED_N >
215 LVARRAY_HOST_DEVICE inline constexpr
216 void checkSizes( T const ( &src )[ INFERRED_M ][ INFERRED_N ] )
217 {
218  static_assert( PROVIDED_M == INFERRED_M, "Expected the first dimension of size PROVIDED_M, got an array of size INFERRED_M." );
219  static_assert( PROVIDED_N == INFERRED_N, "Expected the second dimension of size PROVIDED_N, got an array of size INFERRED_N." );
220  LVARRAY_UNUSED_VARIABLE( src );
221 }
222 
230 template< std::ptrdiff_t ISIZE, typename ARRAY >
231 LVARRAY_HOST_DEVICE inline CONSTEXPR_WITHOUT_BOUNDS_CHECK
232 std::enable_if_t< HasStaticMember_NDIM< ARRAY > >
233 checkSizes( ARRAY const & array )
234 {
235  static_assert( ARRAY::NDIM == 1, "Must be a 1D array." );
236 
237 #ifdef LVARRAY_BOUNDS_CHECK
238  LVARRAY_ERROR_IF_NE( array.size( 0 ), ISIZE );
239 #else
240  LVARRAY_UNUSED_VARIABLE( array );
241 #endif
242 }
243 
252 template< std::ptrdiff_t ISIZE, std::ptrdiff_t JSIZE, typename ARRAY >
253 LVARRAY_HOST_DEVICE inline CONSTEXPR_WITHOUT_BOUNDS_CHECK
254 std::enable_if_t< HasStaticMember_NDIM< ARRAY > >
255 checkSizes( ARRAY const & array )
256 {
257  static_assert( ARRAY::NDIM == 2, "Must be a 1D array." );
258 #ifdef LVARRAY_BOUNDS_CHECK
259  LVARRAY_ERROR_IF_NE( array.size( 0 ), ISIZE );
260  LVARRAY_ERROR_IF_NE( array.size( 1 ), JSIZE );
261 #else
262  LVARRAY_UNUSED_VARIABLE( array );
263 #endif
264 }
265 
266 } // namespace internal
267 
269 template< std::ptrdiff_t ISIZE >
270 constexpr std::ptrdiff_t SYM_SIZE = ( ISIZE * ( ISIZE + 1 ) ) / 2;
271 
278 
286 template< std::ptrdiff_t ISIZE, typename VECTOR >
287 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
288 auto maxAbsoluteEntry( VECTOR && vector )
289 {
290  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
291  internal::checkSizes< ISIZE >( vector );
292 
293  auto maxVal = math::abs( vector[ 0 ] );
294  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
295  {
296  maxVal = math::max( maxVal, math::abs( vector[ i ] ) );
297  }
298 
299  return maxVal;
300 }
301 
310 template< std::ptrdiff_t ISIZE, typename VECTOR >
311 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
312 void fill( VECTOR && vector, std::remove_reference_t< decltype( vector[ 0 ] ) > const value )
313 {
314  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
315  internal::checkSizes< ISIZE >( vector );
316 
317  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
318  {
319  vector[ i ] = value;
320  }
321 }
322 
332 template< std::ptrdiff_t ISIZE, std::ptrdiff_t JSIZE, typename MATRIX >
333 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
334 void fill( MATRIX && matrix, std::remove_reference_t< decltype( matrix[ 0 ][ 0 ] ) > const value )
335 {
336  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
337  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
338  internal::checkSizes< ISIZE, JSIZE >( matrix );
339 
340  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
341  {
342  for( std::ptrdiff_t j = 0; j < JSIZE; ++j )
343  {
344  matrix[ i ][ j ] = value;
345  }
346  }
347 }
348 
358 template< std::ptrdiff_t ISIZE, typename DST_VECTOR, typename SRC_VECTOR >
359 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
360 void copy( DST_VECTOR && LVARRAY_RESTRICT_REF dstVector,
361  SRC_VECTOR const & LVARRAY_RESTRICT_REF srcVector )
362 {
363  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
364  internal::checkSizes< ISIZE >( dstVector );
365  internal::checkSizes< ISIZE >( srcVector );
366 
367  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
368  {
369  dstVector[ i ] = srcVector[ i ];
370  }
371 }
372 
383 template< std::ptrdiff_t ISIZE, std::ptrdiff_t JSIZE, typename DST_MATRIX, typename SRC_MATRIX >
384 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
385 void copy( DST_MATRIX && LVARRAY_RESTRICT_REF dstMatrix,
386  SRC_MATRIX const & LVARRAY_RESTRICT_REF srcMatrix )
387 {
388  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
389  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
390  internal::checkSizes< ISIZE, JSIZE >( dstMatrix );
391  internal::checkSizes< ISIZE, JSIZE >( srcMatrix );
392 
393  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
394  {
395  for( std::ptrdiff_t j = 0; j < JSIZE; ++j )
396  {
397  dstMatrix[ i ][ j ] = srcMatrix[ i ][ j ];
398  }
399  }
400 }
401 
410 template< std::ptrdiff_t ISIZE, typename VECTOR >
411 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
412 void scale( VECTOR && vector, std::remove_reference_t< decltype( vector[ 0 ] ) > const scale )
413 {
414  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
415  internal::checkSizes< ISIZE >( vector );
416 
417  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
418  {
419  vector[ i ] *= scale;
420  }
421 }
422 
432 template< std::ptrdiff_t ISIZE, std::ptrdiff_t JSIZE, typename MATRIX >
433 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
434 void scale( MATRIX && matrix, std::remove_reference_t< decltype( matrix[ 0 ][ 0 ] ) > const scale )
435 {
436  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
437  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
438  internal::checkSizes< ISIZE, JSIZE >( matrix );
439 
440  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
441  {
442  for( std::ptrdiff_t j = 0; j < JSIZE; ++j )
443  {
444  matrix[ i ][ j ] *= scale;
445  }
446  }
447 }
448 
459 template< std::ptrdiff_t ISIZE, typename DST_VECTOR, typename SRC_VECTOR >
460 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
461 void scaledCopy( DST_VECTOR && LVARRAY_RESTRICT_REF dstVector,
462  SRC_VECTOR const & LVARRAY_RESTRICT_REF srcVector,
463  std::remove_reference_t< decltype( srcVector[ 0 ] ) > const scale )
464 {
465  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
466  internal::checkSizes< ISIZE >( dstVector );
467  internal::checkSizes< ISIZE >( srcVector );
468 
469  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
470  {
471  dstVector[ i ] = scale * srcVector[ i ];
472  }
473 }
474 
486 template< std::ptrdiff_t ISIZE, std::ptrdiff_t JSIZE, typename DST_MATRIX, typename SRC_MATRIX >
487 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
488 void scaledCopy( DST_MATRIX && LVARRAY_RESTRICT_REF dstMatrix,
489  SRC_MATRIX const & LVARRAY_RESTRICT_REF srcMatrix,
490  std::remove_reference_t< decltype( srcMatrix[ 0 ][ 0 ] ) > const scale )
491 {
492  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
493  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
494  internal::checkSizes< ISIZE, JSIZE >( dstMatrix );
495  internal::checkSizes< ISIZE, JSIZE >( srcMatrix );
496 
497  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
498  {
499  for( std::ptrdiff_t j = 0; j < JSIZE; ++j )
500  {
501  dstMatrix[ i ][ j ] = scale * srcMatrix[ i ][ j ];
502  }
503  }
504 }
505 
514 template< std::ptrdiff_t M, typename DST_VECTOR >
515 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
516 void addScalar( DST_VECTOR && LVARRAY_RESTRICT_REF dstVector,
517  std::remove_reference_t< decltype( dstVector[0] ) > const value )
518 {
519  static_assert( M > 0, "M must be greater than zero." );
520  internal::checkSizes< M >( dstVector );
521 
522  for( std::ptrdiff_t i = 0; i < M; ++i )
523  {
524  dstVector[ i ] += value;
525  }
526 }
527 
537 template< std::ptrdiff_t ISIZE, typename DST_VECTOR, typename SRC_VECTOR >
538 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
539 void add( DST_VECTOR && LVARRAY_RESTRICT_REF dstVector,
540  SRC_VECTOR const & LVARRAY_RESTRICT_REF srcVector )
541 {
542  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
543  internal::checkSizes< ISIZE >( dstVector );
544  internal::checkSizes< ISIZE >( srcVector );
545 
546  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
547  {
548  dstVector[ i ] += srcVector[ i ];
549  }
550 }
551 
562 template< std::ptrdiff_t ISIZE, std::ptrdiff_t JSIZE, typename DST_MATRIX, typename SRC_MATRIX >
563 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
564 void add( DST_MATRIX && LVARRAY_RESTRICT_REF dstMatrix,
565  SRC_MATRIX const & LVARRAY_RESTRICT_REF srcMatrix )
566 {
567  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
568  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
569  internal::checkSizes< ISIZE, JSIZE >( dstMatrix );
570  internal::checkSizes< ISIZE, JSIZE >( srcMatrix );
571 
572  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
573  {
574  for( std::ptrdiff_t j = 0; j < JSIZE; ++j )
575  {
576  dstMatrix[ i ][ j ] += srcMatrix[ i ][ j ];
577  }
578  }
579 }
580 
590 template< std::ptrdiff_t ISIZE, typename DST_VECTOR, typename SRC_VECTOR >
591 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
592 void subtract( DST_VECTOR && LVARRAY_RESTRICT_REF dstVector,
593  SRC_VECTOR const & LVARRAY_RESTRICT_REF srcVector )
594 {
595  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
596  internal::checkSizes< ISIZE >( dstVector );
597  internal::checkSizes< ISIZE >( srcVector );
598 
599  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
600  {
601  dstVector[ i ] -= srcVector[ i ];
602  }
603 }
604 
615 template< std::ptrdiff_t ISIZE, typename DST_VECTOR, typename SRC_VECTOR >
616 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
617 void scaledAdd( DST_VECTOR && LVARRAY_RESTRICT_REF dstVector,
618  SRC_VECTOR const & LVARRAY_RESTRICT_REF srcVector,
619  std::remove_reference_t< decltype( srcVector[ 0 ] ) > const scale )
620 {
621  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
622  internal::checkSizes< ISIZE >( dstVector );
623  internal::checkSizes< ISIZE >( srcVector );
624 
625  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
626  {
627  dstVector[ i ] = dstVector[ i ] + scale * srcVector[ i ];
628  }
629 }
630 
631 
643 template< std::ptrdiff_t ISIZE, std::ptrdiff_t JSIZE, typename DST_MATRIX, typename SRC_MATRIX >
644 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
645 void scaledAdd( DST_MATRIX && LVARRAY_RESTRICT_REF dstMatrix,
646  SRC_MATRIX const & LVARRAY_RESTRICT_REF srcMatrix,
647  std::remove_reference_t< decltype( srcMatrix[ 0 ][ 0 ] ) > const scale )
648 {
649  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
650  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
651  internal::checkSizes< ISIZE, JSIZE >( dstMatrix );
652  internal::checkSizes< ISIZE, JSIZE >( srcMatrix );
653 
654  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
655  {
656  for( std::ptrdiff_t j = 0; j < JSIZE; ++j )
657  {
658  dstMatrix[ i ][ j ] += scale * srcMatrix[ i ][ j ];
659  }
660  }
661 }
662 
673 template< std::ptrdiff_t ISIZE, typename DST_VECTOR, typename VECTOR_A, typename VECTOR_B >
674 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
675 void hadamardProduct( DST_VECTOR && LVARRAY_RESTRICT_REF dstVector,
676  VECTOR_A const & LVARRAY_RESTRICT_REF vectorA,
677  VECTOR_B const & LVARRAY_RESTRICT_REF vectorB )
678 {
679  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
680  internal::checkSizes< ISIZE >( dstVector );
681  internal::checkSizes< ISIZE >( vectorA );
682  internal::checkSizes< ISIZE >( vectorB );
683 
684  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
685  {
686  dstVector[ i ] = vectorA[ i ] * vectorB[ i ];
687  }
688 }
689 
691 
696 
704 template< std::ptrdiff_t ISIZE, typename VECTOR >
705 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
706 auto l2NormSquared( VECTOR const & vector )
707 {
708  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
709  internal::checkSizes< ISIZE >( vector );
710 
711  auto norm = vector[ 0 ] * vector[ 0 ];
712  for( std::ptrdiff_t i = 1; i < ISIZE; ++i )
713  {
714  norm = norm + vector[ i ] * vector[ i ];
715  }
716  return norm;
717 }
718 
725 template< std::ptrdiff_t ISIZE, typename VECTOR >
726 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
727 auto l2Norm( VECTOR const & vector )
728 {
729  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
730  internal::checkSizes< ISIZE >( vector );
731 
732  return math::sqrt( l2NormSquared< ISIZE >( vector ) );
733 }
734 
742 template< std::ptrdiff_t ISIZE, typename VECTOR >
743 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
744 auto normalize( VECTOR && vector )
745 {
746  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
747  internal::checkSizes< ISIZE >( vector );
748 
749  auto const normInv = math::invSqrt( l2NormSquared< ISIZE >( vector ) );
750  scale< ISIZE >( vector, normInv );
751 
752  return 1 / normInv;
753 }
754 
763 template< std::ptrdiff_t JSIZE, typename VECTOR_A, typename VECTOR_B >
764 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
765 auto AiBi( VECTOR_A const & LVARRAY_RESTRICT_REF vectorA,
766  VECTOR_B const & LVARRAY_RESTRICT_REF vectorB )
767 {
768  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
769  internal::checkSizes< JSIZE >( vectorA );
770  internal::checkSizes< JSIZE >( vectorB );
771 
772  auto result = vectorA[ 0 ] * vectorB[ 0 ];
773  for( std::ptrdiff_t i = 1; i < JSIZE; ++i )
774  {
775  result = result + vectorA[ i ] * vectorB[ i ];
776  }
777 
778  return result;
779 }
780 
790 template< typename DST_VECTOR, typename VECTOR_A, typename VECTOR_B >
791 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
792 void crossProduct( DST_VECTOR && LVARRAY_RESTRICT_REF dstVector,
793  VECTOR_A const & LVARRAY_RESTRICT_REF vectorA,
794  VECTOR_B const & LVARRAY_RESTRICT_REF vectorB )
795 {
796  internal::checkSizes< 3 >( dstVector );
797  internal::checkSizes< 3 >( vectorA );
798  internal::checkSizes< 3 >( vectorB );
799 
800  dstVector[ 0 ] = vectorA[ 1 ] * vectorB[ 2 ] - vectorA[ 2 ] * vectorB[ 1 ];
801  dstVector[ 1 ] = vectorA[ 2 ] * vectorB[ 0 ] - vectorA[ 0 ] * vectorB[ 2 ];
802  dstVector[ 2 ] = vectorA[ 0 ] * vectorB[ 1 ] - vectorA[ 1 ] * vectorB[ 0 ];
803 }
804 
806 
811 
825 template< std::ptrdiff_t ISIZE, std::ptrdiff_t JSIZE, typename DST_MATRIX, typename VECTOR_A, typename VECTOR_B >
826 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
827 void Rij_eq_AiBj( DST_MATRIX && LVARRAY_RESTRICT_REF dstMatrix,
828  VECTOR_A const & LVARRAY_RESTRICT_REF vectorA,
829  VECTOR_B const & LVARRAY_RESTRICT_REF vectorB )
830 {
831  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
832  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
833  internal::checkSizes< ISIZE, JSIZE >( dstMatrix );
834  internal::checkSizes< ISIZE >( vectorA );
835  internal::checkSizes< JSIZE >( vectorB );
836 
837  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
838  {
839  for( std::ptrdiff_t j = 0; j < JSIZE; ++j )
840  {
841  dstMatrix[ i ][ j ] = vectorA[ i ] * vectorB[ j ];
842  }
843  }
844 }
845 
858 template< std::ptrdiff_t JSIZE, std::ptrdiff_t ISIZE, typename DST_MATRIX, typename VECTOR_A, typename VECTOR_B >
859 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
860 void Rij_add_AiBj( DST_MATRIX && LVARRAY_RESTRICT_REF dstMatrix,
861  VECTOR_A const & LVARRAY_RESTRICT_REF vectorA,
862  VECTOR_B const & LVARRAY_RESTRICT_REF vectorB )
863 {
864  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
865  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
866  internal::checkSizes< JSIZE, ISIZE >( dstMatrix );
867  internal::checkSizes< JSIZE >( vectorA );
868  internal::checkSizes< ISIZE >( vectorB );
869 
870  for( std::ptrdiff_t i = 0; i < JSIZE; ++i )
871  {
872  for( std::ptrdiff_t j = 0; j < ISIZE; ++j )
873  {
874  dstMatrix[ i ][ j ] = dstMatrix[ i ][ j ] + vectorA[ i ] * vectorB[ j ];
875  }
876  }
877 }
878 
891 template< std::ptrdiff_t ISIZE, std::ptrdiff_t JSIZE, typename DST_VECTOR, typename MATRIX_A, typename VECTOR_B >
892 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
893 void Ri_eq_AijBj( DST_VECTOR && LVARRAY_RESTRICT_REF dstVector,
894  MATRIX_A const & LVARRAY_RESTRICT_REF matrixA,
895  VECTOR_B const & LVARRAY_RESTRICT_REF vectorB )
896 {
897  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
898  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
899  internal::checkSizes< ISIZE >( dstVector );
900  internal::checkSizes< ISIZE, JSIZE >( matrixA );
901  internal::checkSizes< JSIZE >( vectorB );
902 
903  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
904  {
905  dstVector[ i ] = matrixA[ i ][ 0 ] * vectorB[ 0 ];
906  for( std::ptrdiff_t j = 1; j < JSIZE; ++j )
907  {
908  dstVector[ i ] = dstVector[ i ] + matrixA[ i ][ j ] * vectorB[ j ];
909  }
910  }
911 }
912 
925 template< std::ptrdiff_t ISIZE, std::ptrdiff_t JSIZE, typename DST_VECTOR, typename MATRIX_A, typename VECTOR_B >
926 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
927 void Ri_add_AijBj( DST_VECTOR && LVARRAY_RESTRICT_REF dstVector,
928  MATRIX_A const & LVARRAY_RESTRICT_REF matrixA,
929  VECTOR_B const & LVARRAY_RESTRICT_REF vectorB )
930 {
931  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
932  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
933  internal::checkSizes< ISIZE >( dstVector );
934  internal::checkSizes< ISIZE, JSIZE >( matrixA );
935  internal::checkSizes< JSIZE >( vectorB );
936 
937  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
938  {
939  for( std::ptrdiff_t j = 0; j < JSIZE; ++j )
940  {
941  dstVector[ i ] = dstVector[ i ] + matrixA[ i ][ j ] * vectorB[ j ];
942  }
943  }
944 }
945 
959 template< std::ptrdiff_t ISIZE, std::ptrdiff_t JSIZE, typename DST_VECTOR, typename MATRIX_A, typename VECTOR_B >
960 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
961 void Ri_eq_AjiBj( DST_VECTOR && LVARRAY_RESTRICT_REF dstVector,
962  MATRIX_A const & LVARRAY_RESTRICT_REF matrixA,
963  VECTOR_B const & LVARRAY_RESTRICT_REF vectorB )
964 {
965  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
966  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
967  internal::checkSizes< ISIZE >( dstVector );
968  internal::checkSizes< JSIZE, ISIZE >( matrixA );
969  internal::checkSizes< JSIZE >( vectorB );
970 
971  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
972  {
973  dstVector[ i ] = matrixA[ 0 ][ i ] * vectorB[ 0 ];
974  for( std::ptrdiff_t j = 1; j < JSIZE; ++j )
975  {
976  dstVector[ i ] = dstVector[ i ] + matrixA[ j ][ i ] * vectorB[ j ];
977  }
978  }
979 }
980 
994 template< std::ptrdiff_t ISIZE, std::ptrdiff_t JSIZE, typename DST_VECTOR, typename MATRIX_A, typename VECTOR_B >
995 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
996 void Ri_add_AjiBj( DST_VECTOR && LVARRAY_RESTRICT_REF dstVector,
997  MATRIX_A const & LVARRAY_RESTRICT_REF matrixA,
998  VECTOR_B const & LVARRAY_RESTRICT_REF vectorB )
999 {
1000  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
1001  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
1002  internal::checkSizes< ISIZE >( dstVector );
1003  internal::checkSizes< JSIZE, ISIZE >( matrixA );
1004  internal::checkSizes< JSIZE >( vectorB );
1005 
1006  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
1007  {
1008  for( std::ptrdiff_t j = 0; j < JSIZE; ++j )
1009  {
1010  dstVector[ i ] = dstVector[ i ] + matrixA[ j ][ i ] * vectorB[ j ];
1011  }
1012  }
1013 }
1014 
1016 
1021 
1032 template< std::ptrdiff_t ISIZE, std::ptrdiff_t JSIZE, typename DST_MATRIX, typename SRC_MATRIX >
1033 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
1034 void transpose( DST_MATRIX && LVARRAY_RESTRICT_REF dstMatrix,
1035  SRC_MATRIX const & LVARRAY_RESTRICT_REF srcMatrix )
1036 {
1037  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
1038  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
1039  internal::checkSizes< ISIZE, JSIZE >( dstMatrix );
1040  internal::checkSizes< JSIZE, ISIZE >( srcMatrix );
1041 
1042  for( std::ptrdiff_t i = 0; i < JSIZE; ++i )
1043  {
1044  for( std::ptrdiff_t j = 0; j < ISIZE; ++j )
1045  {
1046  dstMatrix[ j ][ i ] = srcMatrix[ i ][ j ];
1047  }
1048  }
1049 }
1050 
1065 template< std::ptrdiff_t ISIZE,
1066  std::ptrdiff_t JSIZE,
1067  std::ptrdiff_t KSIZE,
1068  typename DST_MATRIX,
1069  typename MATRIX_A,
1070  typename MATRIX_B >
1071 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
1072 void Rij_eq_AikBkj( DST_MATRIX && LVARRAY_RESTRICT_REF dstMatrix,
1073  MATRIX_A const & LVARRAY_RESTRICT_REF matrixA,
1074  MATRIX_B const & LVARRAY_RESTRICT_REF matrixB )
1075 {
1076  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
1077  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
1078  static_assert( KSIZE > 0, "KSIZE must be greater than zero." );
1079  internal::checkSizes< ISIZE, JSIZE >( dstMatrix );
1080  internal::checkSizes< ISIZE, KSIZE >( matrixA );
1081  internal::checkSizes< KSIZE, JSIZE >( matrixB );
1082 
1083 
1084  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
1085  {
1086  for( std::ptrdiff_t j = 0; j < JSIZE; ++j )
1087  {
1088  dstMatrix[ i ][ j ] = matrixA[ i ][ 0 ] * matrixB[ 0 ][ j ];
1089  for( std::ptrdiff_t k = 1; k < KSIZE; ++k )
1090  {
1091  dstMatrix[ i ][ j ] = dstMatrix[ i ][ j ] + matrixA[ i ][ k ] * matrixB[ k ][ j ];
1092  }
1093  }
1094  }
1095 }
1096 
1111 template< std::ptrdiff_t ISIZE,
1112  std::ptrdiff_t JSIZE,
1113  std::ptrdiff_t KSIZE,
1114  typename DST_MATRIX,
1115  typename MATRIX_A,
1116  typename MATRIX_B >
1117 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
1118 void Rij_add_AikBkj( DST_MATRIX && LVARRAY_RESTRICT_REF dstMatrix,
1119  MATRIX_A const & LVARRAY_RESTRICT_REF matrixA,
1120  MATRIX_B const & LVARRAY_RESTRICT_REF matrixB )
1121 {
1122  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
1123  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
1124  static_assert( KSIZE > 0, "KSIZE must be greater than zero." );
1125  internal::checkSizes< ISIZE, JSIZE >( dstMatrix );
1126  internal::checkSizes< ISIZE, KSIZE >( matrixA );
1127  internal::checkSizes< KSIZE, JSIZE >( matrixB );
1128 
1129 
1130  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
1131  {
1132  for( std::ptrdiff_t j = 0; j < JSIZE; ++j )
1133  {
1134  for( std::ptrdiff_t k = 0; k < KSIZE; ++k )
1135  {
1136  dstMatrix[ i ][ j ] = dstMatrix[ i ][ j ] + matrixA[ i ][ k ] * matrixB[ k ][ j ];
1137  }
1138  }
1139  }
1140 }
1141 
1156 template< std::ptrdiff_t ISIZE,
1157  std::ptrdiff_t JSIZE,
1158  std::ptrdiff_t KSIZE,
1159  typename DST_MATRIX,
1160  typename MATRIX_A,
1161  typename MATRIX_B >
1162 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
1163 void Rij_eq_AikBjk( DST_MATRIX && LVARRAY_RESTRICT_REF dstMatrix,
1164  MATRIX_A const & LVARRAY_RESTRICT_REF matrixA,
1165  MATRIX_B const & LVARRAY_RESTRICT_REF matrixB )
1166 {
1167  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
1168  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
1169  static_assert( KSIZE > 0, "KSIZE must be greater than zero." );
1170  internal::checkSizes< ISIZE, JSIZE >( dstMatrix );
1171  internal::checkSizes< ISIZE, KSIZE >( matrixA );
1172  internal::checkSizes< JSIZE, KSIZE >( matrixB );
1173 
1174 
1175  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
1176  {
1177  for( std::ptrdiff_t j = 0; j < JSIZE; ++j )
1178  {
1179  dstMatrix[ i ][ j ] = matrixA[ i ][ 0 ] * matrixB[ j ][ 0 ];
1180  for( std::ptrdiff_t k = 1; k < KSIZE; ++k )
1181  {
1182  dstMatrix[ i ][ j ] = dstMatrix[ i ][ j ] + matrixA[ i ][ k ] * matrixB[ j ][ k ];
1183  }
1184  }
1185  }
1186 }
1187 
1202 template< std::ptrdiff_t ISIZE,
1203  std::ptrdiff_t JSIZE,
1204  std::ptrdiff_t KSIZE,
1205  typename DST_MATRIX,
1206  typename MATRIX_A,
1207  typename MATRIX_B >
1208 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
1209 void Rij_add_AikBjk( DST_MATRIX && LVARRAY_RESTRICT_REF dstMatrix,
1210  MATRIX_A const & LVARRAY_RESTRICT_REF matrixA,
1211  MATRIX_B const & LVARRAY_RESTRICT_REF matrixB )
1212 {
1213  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
1214  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
1215  static_assert( KSIZE > 0, "KSIZE must be greater than zero." );
1216  internal::checkSizes< ISIZE, JSIZE >( dstMatrix );
1217  internal::checkSizes< ISIZE, KSIZE >( matrixA );
1218  internal::checkSizes< JSIZE, KSIZE >( matrixB );
1219 
1220  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
1221  {
1222  for( std::ptrdiff_t j = 0; j < JSIZE; ++j )
1223  {
1224  for( std::ptrdiff_t k = 0; k < KSIZE; ++k )
1225  {
1226  dstMatrix[ i ][ j ] = dstMatrix[ i ][ j ] + matrixA[ i ][ k ] * matrixB[ j ][ k ];
1227  }
1228  }
1229  }
1230 }
1231 
1243 template< std::ptrdiff_t ISIZE,
1244  std::ptrdiff_t JSIZE,
1245  typename DST_MATRIX,
1246  typename MATRIX_A >
1247 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
1248 void Rij_add_AikAjk( DST_MATRIX && LVARRAY_RESTRICT_REF dstMatrix,
1249  MATRIX_A const & LVARRAY_RESTRICT_REF matrixA )
1250 {
1251  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
1252  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
1253  internal::checkSizes< ISIZE, ISIZE >( dstMatrix );
1254  internal::checkSizes< ISIZE, JSIZE >( matrixA );
1255 
1256  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
1257  {
1258  for( std::ptrdiff_t j = 0; j < ISIZE; ++j )
1259  {
1260  for( std::ptrdiff_t k = 0; k < JSIZE; ++k )
1261  {
1262  dstMatrix[ i ][ j ] = dstMatrix[ i ][ j ] + matrixA[ i ][ k ] * matrixA[ j ][ k ];
1263  }
1264  }
1265  }
1266 }
1267 
1282 template< std::ptrdiff_t ISIZE,
1283  std::ptrdiff_t JSIZE,
1284  std::ptrdiff_t KSIZE,
1285  typename DST_MATRIX,
1286  typename MATRIX_A,
1287  typename MATRIX_B >
1288 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
1289 void Rij_eq_AkiBkj( DST_MATRIX && LVARRAY_RESTRICT_REF dstMatrix,
1290  MATRIX_A const & LVARRAY_RESTRICT_REF matrixA,
1291  MATRIX_B const & LVARRAY_RESTRICT_REF matrixB )
1292 {
1293  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
1294  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
1295  static_assert( KSIZE > 0, "KSIZE must be greater than zero." );
1296  internal::checkSizes< ISIZE, JSIZE >( dstMatrix );
1297  internal::checkSizes< KSIZE, ISIZE >( matrixA );
1298  internal::checkSizes< KSIZE, JSIZE >( matrixB );
1299 
1300  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
1301  {
1302  for( std::ptrdiff_t j = 0; j < JSIZE; ++j )
1303  {
1304  dstMatrix[ i ][ j ] = matrixA[ 0 ][ i ] * matrixB[ 0 ][ j ];
1305  for( std::ptrdiff_t k = 1; k < KSIZE; ++k )
1306  {
1307  dstMatrix[ i ][ j ] = dstMatrix[ i ][ j ] + matrixA[ k ][ i ] * matrixB[ k ][ j ];
1308  }
1309  }
1310  }
1311 }
1312 
1324 template< std::ptrdiff_t ISIZE,
1325  std::ptrdiff_t JSIZE,
1326  typename DST_MATRIX,
1327  typename MATRIX_A >
1328 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
1329 void Rij_eq_AkiAkj( DST_MATRIX && LVARRAY_RESTRICT_REF dstMatrix,
1330  MATRIX_A const & LVARRAY_RESTRICT_REF matrixA )
1331 {
1332  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
1333  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
1334  internal::checkSizes< ISIZE, ISIZE >( dstMatrix );
1335  internal::checkSizes< JSIZE, ISIZE >( matrixA );
1336 
1337  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
1338  {
1339  for( std::ptrdiff_t j = 0; j < ISIZE; ++j )
1340  {
1341  dstMatrix[ i ][ j ] = matrixA[ 0 ][ i ] * matrixA[ 0 ][ j ];
1342  for( std::ptrdiff_t k = 1; k < JSIZE; ++k )
1343  {
1344  dstMatrix[ i ][ j ] = dstMatrix[ i ][ j ] + matrixA[ k ][ i ] * matrixA[ k ][ j ];
1345  }
1346  }
1347  }
1348 }
1349 
1364 template< std::ptrdiff_t ISIZE,
1365  std::ptrdiff_t JSIZE,
1366  std::ptrdiff_t KSIZE,
1367  typename DST_MATRIX,
1368  typename MATRIX_A,
1369  typename MATRIX_B >
1370 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
1371 void Rij_add_AkiBkj( DST_MATRIX && LVARRAY_RESTRICT_REF dstMatrix,
1372  MATRIX_A const & LVARRAY_RESTRICT_REF matrixA,
1373  MATRIX_B const & LVARRAY_RESTRICT_REF matrixB )
1374 {
1375  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
1376  static_assert( JSIZE > 0, "JSIZE must be greater than zero." );
1377  static_assert( KSIZE > 0, "KSIZE must be greater than zero." );
1378  internal::checkSizes< ISIZE, JSIZE >( dstMatrix );
1379  internal::checkSizes< KSIZE, ISIZE >( matrixA );
1380  internal::checkSizes< KSIZE, JSIZE >( matrixB );
1381 
1382  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
1383  {
1384  for( std::ptrdiff_t j = 0; j < JSIZE; ++j )
1385  {
1386  for( std::ptrdiff_t k = 0; k < KSIZE; ++k )
1387  {
1388  dstMatrix[ i ][ j ] = dstMatrix[ i ][ j ] + matrixA[ k ][ i ] * matrixB[ k ][ j ];
1389  }
1390  }
1391  }
1392 }
1393 
1395 
1400 
1408 template< std::ptrdiff_t ISIZE, typename MATRIX >
1409 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
1410 void transpose( MATRIX && LVARRAY_RESTRICT_REF matrix )
1411 {
1412  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
1413  internal::checkSizes< ISIZE, ISIZE >( matrix );
1414 
1415  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
1416  {
1417  for( std::ptrdiff_t j = i + 1; j < ISIZE; ++j )
1418  {
1419  auto const entryIJ = matrix[ i ][ j ];
1420  matrix[ i ][ j ] = matrix[ j ][ i ];
1421  matrix[ j ][ i ] = entryIJ;
1422  }
1423  }
1424 }
1425 
1434 template< std::ptrdiff_t ISIZE, typename MATRIX >
1435 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
1436 void addIdentity( MATRIX && matrix, std::remove_reference_t< decltype( matrix[ 0 ][ 0 ] ) > const scale )
1437 {
1438  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
1439  internal::checkSizes< ISIZE, ISIZE >( matrix );
1440 
1441  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
1442  {
1443  matrix[ i ][ i ] += scale;
1444  }
1445 }
1446 
1453 template< std::ptrdiff_t ISIZE, typename MATRIX >
1454 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
1455 auto trace( MATRIX const & matrix )
1456 {
1457  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
1458  internal::checkSizes< ISIZE, ISIZE >( matrix );
1459 
1460  auto trace = matrix[ 0 ][ 0 ];
1461  for( std::ptrdiff_t i = 1; i < ISIZE; ++i )
1462  {
1463  trace += matrix[ i ][ i ];
1464  }
1465 
1466  return trace;
1467 }
1468 
1470 
1475 
1485 template< std::ptrdiff_t ISIZE, typename SYM_MATRIX >
1486 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
1487 void symAddIdentity( SYM_MATRIX && symMatrix, std::remove_reference_t< decltype( symMatrix[ 0 ] ) > const scale )
1488 {
1489  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
1490  internal::checkSizes< SYM_SIZE< ISIZE > >( symMatrix );
1491 
1492  for( std::ptrdiff_t i = 0; i < ISIZE; ++i )
1493  {
1494  symMatrix[ i ] += scale;
1495  }
1496 }
1497 
1504 template< std::ptrdiff_t ISIZE, typename SYM_MATRIX >
1505 LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK inline
1506 auto symTrace( SYM_MATRIX const & symMatrix )
1507 {
1508  static_assert( ISIZE > 0, "ISIZE must be greater than zero." );
1509  internal::checkSizes< SYM_SIZE< ISIZE > >( symMatrix );
1510 
1511  auto trace = symMatrix[ 0 ];
1512  for( std::ptrdiff_t i = 1; i < ISIZE; ++i )
1513  {
1514  trace += symMatrix[ i ];
1515  }
1516 
1517  return trace;
1518 }
1519 
1521 
1522 } // namespace tensorOps
1523 } // namespace LvArray
LVARRAY_UNUSED_VARIABLE
#define LVARRAY_UNUSED_VARIABLE(X)
Mark X as an unused variable, used to silence compiler warnings.
Definition: Macros.hpp:51
LvArray::tensorOps::Rij_add_AikBkj
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void Rij_add_AikBkj(DST_MATRIX &&LVARRAY_RESTRICT_REF dstMatrix, MATRIX_A const &LVARRAY_RESTRICT_REF matrixA, MATRIX_B const &LVARRAY_RESTRICT_REF matrixB)
Multiply matrixA with matrixB and add the result to dstMatrix.
Definition: genericTensorOps.hpp:1118
LvArray::tensorOps::subtract
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void subtract(DST_VECTOR &&LVARRAY_RESTRICT_REF dstVector, SRC_VECTOR const &LVARRAY_RESTRICT_REF srcVector)
Subtract srcVector from dstVector.
Definition: genericTensorOps.hpp:592
LvArray::tensorOps::scaledCopy
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void scaledCopy(DST_VECTOR &&LVARRAY_RESTRICT_REF dstVector, SRC_VECTOR const &LVARRAY_RESTRICT_REF srcVector, std::remove_reference_t< decltype(srcVector[0]) > const scale)
Copy srcVector scaled by scale into dstVector.
Definition: genericTensorOps.hpp:461
LvArray::tensorOps::SYM_SIZE
constexpr std::ptrdiff_t SYM_SIZE
The size of a symmetric MxM matrix in Voigt notation.
Definition: genericTensorOps.hpp:270
LvArray::tensorOps::scaledAdd
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void scaledAdd(DST_VECTOR &&LVARRAY_RESTRICT_REF dstVector, SRC_VECTOR const &LVARRAY_RESTRICT_REF srcVector, std::remove_reference_t< decltype(srcVector[0]) > const scale)
Add srcVector scaled by scale to dstVector.
Definition: genericTensorOps.hpp:617
LvArray::tensorOps::add
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void add(DST_VECTOR &&LVARRAY_RESTRICT_REF dstVector, SRC_VECTOR const &LVARRAY_RESTRICT_REF srcVector)
Add srcVector to dstVector.
Definition: genericTensorOps.hpp:539
LvArray::tensorOps::fill
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void fill(VECTOR &&vector, std::remove_reference_t< decltype(vector[0]) > const value)
Set the entries of vector to value.
Definition: genericTensorOps.hpp:312
HAS_STATIC_MEMBER
#define HAS_STATIC_MEMBER(NAME)
Macro that expands to a static constexpr bool templated on a type that is only true when the type has...
Definition: typeManipulation.hpp:109
LvArray::tensorOps::Rij_add_AikBjk
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void Rij_add_AikBjk(DST_MATRIX &&LVARRAY_RESTRICT_REF dstMatrix, MATRIX_A const &LVARRAY_RESTRICT_REF matrixA, MATRIX_B const &LVARRAY_RESTRICT_REF matrixB)
Multiply matrixA with the transpose of matrixB and put the result into dstMatrix. ...
Definition: genericTensorOps.hpp:1209
LvArray::tensorOps::Rij_eq_AiBj
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void Rij_eq_AiBj(DST_MATRIX &&LVARRAY_RESTRICT_REF dstMatrix, VECTOR_A const &LVARRAY_RESTRICT_REF vectorA, VECTOR_B const &LVARRAY_RESTRICT_REF vectorB)
Perform the outer product of vectorA and vectorB writing the result to dstMatrix. ...
Definition: genericTensorOps.hpp:827
LvArray::tensorOps::Rij_add_AkiBkj
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void Rij_add_AkiBkj(DST_MATRIX &&LVARRAY_RESTRICT_REF dstMatrix, MATRIX_A const &LVARRAY_RESTRICT_REF matrixA, MATRIX_B const &LVARRAY_RESTRICT_REF matrixB)
Multiply the transpose of matrixA with matrixB and add the result into dstMatrix. ...
Definition: genericTensorOps.hpp:1371
LvArray::tensorOps::Ri_add_AjiBj
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void Ri_add_AjiBj(DST_VECTOR &&LVARRAY_RESTRICT_REF dstVector, MATRIX_A const &LVARRAY_RESTRICT_REF matrixA, VECTOR_B const &LVARRAY_RESTRICT_REF vectorB)
Perform the matrix vector multiplication of the transpose of matrixA and vectorB adding the result to...
Definition: genericTensorOps.hpp:996
LvArray::tensorOps::hadamardProduct
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void hadamardProduct(DST_VECTOR &&LVARRAY_RESTRICT_REF dstVector, VECTOR_A const &LVARRAY_RESTRICT_REF vectorA, VECTOR_B const &LVARRAY_RESTRICT_REF vectorB)
Multiply the elements of vectorA and vectorB putting the result into dstVector.
Definition: genericTensorOps.hpp:675
LvArray::tensorOps::trace
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK auto trace(MATRIX const &matrix)
Definition: genericTensorOps.hpp:1455
LvArray::tensorOps::l2Norm
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK auto l2Norm(VECTOR const &vector)
Definition: genericTensorOps.hpp:727
LvArray::tensorOps::copy
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void copy(DST_VECTOR &&LVARRAY_RESTRICT_REF dstVector, SRC_VECTOR const &LVARRAY_RESTRICT_REF srcVector)
Copy srcVector into dstVector.
Definition: genericTensorOps.hpp:360
CONSTEXPR_WITHOUT_BOUNDS_CHECK
#define CONSTEXPR_WITHOUT_BOUNDS_CHECK
Expands to constexpr when array bound checking is disabled.
Definition: Macros.hpp:609
LvArray::tensorOps::Rij_add_AiBj
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void Rij_add_AiBj(DST_MATRIX &&LVARRAY_RESTRICT_REF dstMatrix, VECTOR_A const &LVARRAY_RESTRICT_REF vectorA, VECTOR_B const &LVARRAY_RESTRICT_REF vectorB)
Perform the outer product of vectorA and vectorB adding the result to dstMatrix.
Definition: genericTensorOps.hpp:860
ArraySlice.hpp
Contains the implementation of LvArray::ArraySlice.
LvArray::tensorOps::l2NormSquared
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK auto l2NormSquared(VECTOR const &vector)
Definition: genericTensorOps.hpp:706
LvArray::tensorOps::Ri_eq_AjiBj
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void Ri_eq_AjiBj(DST_VECTOR &&LVARRAY_RESTRICT_REF dstVector, MATRIX_A const &LVARRAY_RESTRICT_REF matrixA, VECTOR_B const &LVARRAY_RESTRICT_REF vectorB)
Perform the matrix vector multiplication of the transpose of matrixA and vectorB writing the result t...
Definition: genericTensorOps.hpp:961
LvArray::tensorOps::transpose
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void transpose(DST_MATRIX &&LVARRAY_RESTRICT_REF dstMatrix, SRC_MATRIX const &LVARRAY_RESTRICT_REF srcMatrix)
Store the transpose of the NxM matrix srcMatrix in dstMatrix.
Definition: genericTensorOps.hpp:1034
LvArray::tensorOps::Rij_eq_AikBjk
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void Rij_eq_AikBjk(DST_MATRIX &&LVARRAY_RESTRICT_REF dstMatrix, MATRIX_A const &LVARRAY_RESTRICT_REF matrixA, MATRIX_B const &LVARRAY_RESTRICT_REF matrixB)
Multiply matrixA with the transpose of matrixB and put the result into dstMatrix. ...
Definition: genericTensorOps.hpp:1163
math.hpp
Contains some portable math functions.
LvArray::tensorOps::maxAbsoluteEntry
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK auto maxAbsoluteEntry(VECTOR &&vector)
Definition: genericTensorOps.hpp:288
LvArray::tensorOps::AiBi
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK auto AiBi(VECTOR_A const &LVARRAY_RESTRICT_REF vectorA, VECTOR_B const &LVARRAY_RESTRICT_REF vectorB)
Definition: genericTensorOps.hpp:765
LvArray::tensorOps::addScalar
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void addScalar(DST_VECTOR &&LVARRAY_RESTRICT_REF dstVector, std::remove_reference_t< decltype(dstVector[0]) > const value)
Add value to dstVector.
Definition: genericTensorOps.hpp:516
LvArray::tensorOps::symAddIdentity
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void symAddIdentity(SYM_MATRIX &&symMatrix, std::remove_reference_t< decltype(symMatrix[0]) > const scale)
Add scale times the identity matrix to symMatrix.
Definition: genericTensorOps.hpp:1487
LvArray::tensorOps::symTrace
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK auto symTrace(SYM_MATRIX const &symMatrix)
Definition: genericTensorOps.hpp:1506
LvArray::tensorOps::Rij_eq_AikBkj
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void Rij_eq_AikBkj(DST_MATRIX &&LVARRAY_RESTRICT_REF dstMatrix, MATRIX_A const &LVARRAY_RESTRICT_REF matrixA, MATRIX_B const &LVARRAY_RESTRICT_REF matrixB)
Multiply matrixA with matrixB and put the result into dstMatrix.
Definition: genericTensorOps.hpp:1072
LvArray::math::sqrt
LVARRAY_HOST_DEVICE float sqrt(float const x)
Definition: math.hpp:461
LvArray
The top level namespace.
Definition: Array.hpp:24
LvArray::tensorOps::Rij_eq_AkiBkj
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void Rij_eq_AkiBkj(DST_MATRIX &&LVARRAY_RESTRICT_REF dstMatrix, MATRIX_A const &LVARRAY_RESTRICT_REF matrixA, MATRIX_B const &LVARRAY_RESTRICT_REF matrixB)
Multiply the transpose of matrixA with matrixB and put the result into dstMatrix. ...
Definition: genericTensorOps.hpp:1289
LvArray::tensorOps::crossProduct
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void crossProduct(DST_VECTOR &&LVARRAY_RESTRICT_REF dstVector, VECTOR_A const &LVARRAY_RESTRICT_REF vectorA, VECTOR_B const &LVARRAY_RESTRICT_REF vectorB)
Compute the cross product of vectorA and vectorB and put it in dstVector.
Definition: genericTensorOps.hpp:792
LvArray::tensorOps::addIdentity
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void addIdentity(MATRIX &&matrix, std::remove_reference_t< decltype(matrix[0][0]) > const scale)
Add scale times the identity matrix to matrix.
Definition: genericTensorOps.hpp:1436
LvArray::math::abs
LVARRAY_HOST_DEVICE constexpr T abs(T const x)
Definition: math.hpp:402
LvArray::tensorOps::Rij_eq_AkiAkj
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void Rij_eq_AkiAkj(DST_MATRIX &&LVARRAY_RESTRICT_REF dstMatrix, MATRIX_A const &LVARRAY_RESTRICT_REF matrixA)
Multiply the transpose of matrixA with matrixA and put the result into dstMatrix. ...
Definition: genericTensorOps.hpp:1329
LvArray::math::max
LVARRAY_HOST_DEVICE constexpr std::enable_if_t< std::is_arithmetic< T >::value, T > max(T const a, T const b)
Definition: math.hpp:311
LvArray::tensorOps::scale
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void scale(VECTOR &&vector, std::remove_reference_t< decltype(vector[0]) > const scale)
Multiply the entries of vector by scale.
Definition: genericTensorOps.hpp:412
LvArray::math::invSqrt
LVARRAY_HOST_DEVICE float invSqrt(float const x)
Definition: math.hpp:503
LvArray::tensorOps::normalize
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK auto normalize(VECTOR &&vector)
Scale vector to a unit vector.
Definition: genericTensorOps.hpp:744
LvArray::tensorOps::Rij_add_AikAjk
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void Rij_add_AikAjk(DST_MATRIX &&LVARRAY_RESTRICT_REF dstMatrix, MATRIX_A const &LVARRAY_RESTRICT_REF matrixA)
Multiply matrixA with the transpose of itself and put the result into dstMatrix.
Definition: genericTensorOps.hpp:1248
LVARRAY_ERROR_IF_NE
#define LVARRAY_ERROR_IF_NE(lhs, rhs)
Raise a hard error if two values are not equal.
Definition: Macros.hpp:321
LvArray::tensorOps::internal::checkSizes
LVARRAY_HOST_DEVICE constexpr std::enable_if_t< HasStaticMember_SIZE< T > > checkSizes(T const &src)
Verify at compile time that the size of a user-provided type is as expected.
Definition: genericTensorOps.hpp:178
LvArray::tensorOps::Ri_add_AijBj
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void Ri_add_AijBj(DST_VECTOR &&LVARRAY_RESTRICT_REF dstVector, MATRIX_A const &LVARRAY_RESTRICT_REF matrixA, VECTOR_B const &LVARRAY_RESTRICT_REF vectorB)
Perform the matrix vector multiplication of matrixA and vectorB adding the result to dstVector...
Definition: genericTensorOps.hpp:927
LVARRAY_HOST_DEVICE
#define LVARRAY_HOST_DEVICE
Mark a function for both host and device usage.
Definition: Macros.hpp:549
LvArray::tensorOps::Ri_eq_AijBj
LVARRAY_HOST_DEVICE CONSTEXPR_WITHOUT_BOUNDS_CHECK void Ri_eq_AijBj(DST_VECTOR &&LVARRAY_RESTRICT_REF dstVector, MATRIX_A const &LVARRAY_RESTRICT_REF matrixA, VECTOR_B const &LVARRAY_RESTRICT_REF vectorB)
Perform the matrix vector multiplication of matrixA and vectorB writing the result to dstVector...
Definition: genericTensorOps.hpp:893
Generated by   doxygen 1.8.13