libMesh
Public Types | Public Member Functions | Protected Member Functions | Protected Attributes | Friends | List of all members
libMesh::TypeTensor< T > Class Template Reference

This class defines a tensor in LIBMESH_DIM dimensional space of type T. More...

#include <tensor_tools.h>

Inheritance diagram for libMesh::TypeTensor< T >:
[legend]

Public Types

typedef T value_type
 Helper typedef for C++98 generic programming. More...
 
typedef std::tuple< unsigned int, unsigned intindex_type
 Helper typedef for generic index programming. More...
 

Public Member Functions

 TypeTensor ()
 Empty constructor. More...
 
template<typename T2 >
 TypeTensor (const TypeTensor< T2 > &p)
 Copy-constructor. More...
 
 ~TypeTensor ()
 Destructor. More...
 
template<typename T2 >
void assign (const TypeTensor< T2 > &)
 Assign to this tensor without creating a temporary. More...
 
template<typename Scalar >
boostcopy::enable_if_c< ScalarTraits< Scalar >::value, TypeTensor & >::type operator= (const Scalar &libmesh_dbg_var(p))
 Assignment-from-scalar operator. More...
 
const T & operator() (const unsigned int i, const unsigned int j) const
 
T & operator() (const unsigned int i, const unsigned int j)
 
ConstTypeTensorColumn< T > slice (const unsigned int i) const
 
TypeTensorColumn< T > slice (const unsigned int i)
 
TypeVector< T > row (const unsigned int r) const
 
TypeVector< T > column (const unsigned int r) const
 
template<typename T2 >
TypeTensor< typename CompareTypes< T, T2 >::supertype > operator+ (const TypeTensor< T2 > &) const
 Add another tensor to this tensor. More...
 
template<typename T2 >
const TypeTensor< T > & operator+= (const TypeTensor< T2 > &)
 Add to this tensor. More...
 
template<typename T2 >
void add (const TypeTensor< T2 > &)
 Add to this tensor without creating a temporary. More...
 
template<typename T2 >
void add_scaled (const TypeTensor< T2 > &, const T &)
 Add a scaled tensor to this tensor without creating a temporary. More...
 
template<typename T2 >
TypeTensor< typename CompareTypes< T, T2 >::supertype > operator- (const TypeTensor< T2 > &) const
 Subtract a tensor from this tensor. More...
 
template<typename T2 >
const TypeTensor< T > & operator-= (const TypeTensor< T2 > &)
 Subtract from this tensor. More...
 
template<typename T2 >
void subtract (const TypeTensor< T2 > &)
 Subtract from this tensor without creating a temporary. More...
 
template<typename T2 >
void subtract_scaled (const TypeTensor< T2 > &, const T &)
 Subtract a scaled value from this tensor without creating a temporary. More...
 
TypeTensor< T > operator- () const
 
template<typename Scalar >
auto operator* (const Scalar &scalar) const -> typename boostcopy::enable_if_c< ScalarTraits< Scalar >::value, TypeTensor< decltype(T() *scalar)>>::type
 Multiply this tensor by a scalar value. More...
 
template<typename Scalar , typename boostcopy::enable_if_c< ScalarTraits< Scalar >::value, int >::type = 0>
const TypeTensor< T > & operator*= (const Scalar &factor)
 Multiply this tensor by a scalar value in place. More...
 
template<typename Scalar >
boostcopy::enable_if_c< ScalarTraits< Scalar >::value, TypeTensor< typename CompareTypes< T, Scalar >::supertype > >::type operator/ (const Scalar &) const
 Divide each entry of this tensor by a scalar value. More...
 
const TypeTensor< T > & operator/= (const T &)
 Divide each entry of this tensor by a scalar value. More...
 
template<typename T2 >
TypeTensor< typename CompareTypes< T, T2 >::supertype > operator* (const TypeTensor< T2 > &) const
 Multiply 2 tensors together, i.e. More...
 
template<typename T2 >
const TypeTensor< T > & operator*= (const TypeTensor< T2 > &)
 Multiply this tensor by a tensor value in place. More...
 
template<typename T2 >
CompareTypes< T, T2 >::supertype contract (const TypeTensor< T2 > &) const
 Multiply 2 tensors together to return a scalar, i.e. More...
 
template<typename T2 >
TypeVector< typename CompareTypes< T, T2 >::supertype > operator* (const TypeVector< T2 > &) const
 Right-multiply this tensor by a vector, i.e. More...
 
template<typename T2 >
TypeVector< typename CompareTypes< T, T2 >::supertype > left_multiply (const TypeVector< T2 > &p) const
 Left-multiply this tensor by a vector, i.e. More...
 
TypeTensor< T > transpose () const
 
TypeTensor< T > inverse () const
 
void solve (const TypeVector< T > &b, TypeVector< T > &x) const
 Solve the 2x2 or 3x3 system of equations A*x = b for x by directly inverting A and multiplying it by b. More...
 
auto norm () const -> decltype(std::norm(T()))
 
auto norm_sq () const -> decltype(std::norm(T()))
 
bool is_zero () const
 
det () const
 
tr () const
 
void zero ()
 Set all entries of the tensor to 0. More...
 
bool operator== (const TypeTensor< T > &rhs) const
 
bool operator< (const TypeTensor< T > &rhs) const
 
bool operator> (const TypeTensor< T > &rhs) const
 
void print (std::ostream &os=libMesh::out) const
 Formatted print to a stream which defaults to libMesh::out. More...
 
void write_unformatted (std::ostream &out_stream, const bool newline=true) const
 Unformatted print to the stream out. More...
 
template<>
bool operator< (const TypeTensor< Real > &rhs) const
 
template<>
bool operator> (const TypeTensor< Real > &rhs) const
 
template<>
bool operator< (const TypeTensor< Complex > &rhs) const
 
template<>
bool operator> (const TypeTensor< Complex > &rhs) const
 

Protected Member Functions

 TypeTensor (const T &xx, const T &xy=0, const T &xz=0, const T &yx=0, const T &yy=0, const T &yz=0, const T &zx=0, const T &zy=0, const T &zz=0)
 Constructor-from-T. More...
 
template<typename Scalar >
 TypeTensor (const Scalar &xx, const Scalar &xy=0, const Scalar &xz=0, const Scalar &yx=0, const Scalar &yy=0, const Scalar &yz=0, const Scalar &zx=0, const Scalar &zy=0, typename boostcopy::enable_if_c< ScalarTraits< Scalar >::value, const Scalar >::type &zz=0)
 Constructor-from-Scalar. More...
 
template<typename T2 >
 TypeTensor (const TypeVector< T2 > &vx)
 Constructor. More...
 
template<typename T2 >
 TypeTensor (const TypeVector< T2 > &vx, const TypeVector< T2 > &vy)
 
template<typename T2 >
 TypeTensor (const TypeVector< T2 > &vx, const TypeVector< T2 > &vy, const TypeVector< T2 > &vz)
 

Protected Attributes

_coords [LIBMESH_DIM *LIBMESH_DIM]
 The coordinates of the TypeTensor. More...
 

Friends

template<typename T2 >
class TypeTensor
 
std::ostream & operator<< (std::ostream &os, const TypeTensor< T > &t)
 Formatted print as above but supports the syntax: More...
 

Detailed Description

template<typename T>
class libMesh::TypeTensor< T >

This class defines a tensor in LIBMESH_DIM dimensional space of type T.

T may either be Real or Complex.

Author
Roy Stogner
Date
2004

Definition at line 36 of file tensor_tools.h.

Member Typedef Documentation

◆ index_type

template<typename T>
typedef std::tuple<unsigned int, unsigned int> libMesh::TypeTensor< T >::index_type

Helper typedef for generic index programming.

Definition at line 125 of file type_tensor.h.

◆ value_type

template<typename T>
typedef T libMesh::TypeTensor< T >::value_type

Helper typedef for C++98 generic programming.

Definition at line 120 of file type_tensor.h.

Constructor & Destructor Documentation

◆ TypeTensor() [1/7]

template<typename T >
libMesh::TypeTensor< T >::TypeTensor ( )
inline

Empty constructor.

Gives the tensor 0 in LIBMESH_DIM dimensions.

Definition at line 510 of file type_tensor.h.

511 {
512  _coords[0] = {};
513 
514 #if LIBMESH_DIM > 1
515  _coords[1] = {};
516  _coords[2] = {};
517  _coords[3] = {};
518 #endif
519 
520 #if LIBMESH_DIM > 2
521  _coords[4] = {};
522  _coords[5] = {};
523  _coords[6] = {};
524  _coords[7] = {};
525  _coords[8] = {};
526 #endif
527 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ TypeTensor() [2/7]

template<typename T >
libMesh::TypeTensor< T >::TypeTensor ( const T &  xx,
const T &  xy = 0,
const T &  xz = 0,
const T &  yx = 0,
const T &  yy = 0,
const T &  yz = 0,
const T &  zx = 0,
const T &  zy = 0,
const T &  zz = 0 
)
inlineexplicitprotected

Constructor-from-T.

By default sets higher dimensional entries to 0. This is a poor constructor for 2D tensors - if the default arguments are to be overridden it requires that the "xz = 0." etc. arguments also be given explicitly.

Definition at line 533 of file type_tensor.h.

References libMesh::libmesh_ignore().

542 {
543  _coords[0] = xx;
544 
545 #if LIBMESH_DIM == 2
546  _coords[1] = xy;
547  _coords[2] = yx;
548  _coords[3] = yy;
549 #elif LIBMESH_DIM == 1
550  libmesh_assert_equal_to (xy, 0);
551  libmesh_assert_equal_to (yx, 0);
552  libmesh_assert_equal_to (yy, 0);
553  libmesh_ignore(xy, yx, yy);
554 #endif
555 
556 #if LIBMESH_DIM == 3
557  _coords[1] = xy;
558  _coords[2] = xz;
559  _coords[3] = yx;
560  _coords[4] = yy;
561  _coords[5] = yz;
562  _coords[6] = zx;
563  _coords[7] = zy;
564  _coords[8] = zz;
565 #else
566  libmesh_assert_equal_to (xz, 0);
567  libmesh_assert_equal_to (yz, 0);
568  libmesh_assert_equal_to (zx, 0);
569  libmesh_assert_equal_to (zy, 0);
570  libmesh_assert_equal_to (zz, 0);
571  libmesh_ignore(xz, yz, zx, zy, zz);
572 #endif
573 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444
void libmesh_ignore(const Args &...)

◆ TypeTensor() [3/7]

template<typename T >
template<typename Scalar >
libMesh::TypeTensor< T >::TypeTensor ( const Scalar &  xx,
const Scalar &  xy = 0,
const Scalar &  xz = 0,
const Scalar &  yx = 0,
const Scalar &  yy = 0,
const Scalar &  yz = 0,
const Scalar &  zx = 0,
const Scalar &  zy = 0,
typename boostcopy::enable_if_c< ScalarTraits< Scalar >::value, const Scalar >::type &  zz = 0 
)
inlineexplicitprotected

Constructor-from-Scalar.

Definition at line 579 of file type_tensor.h.

References libMesh::libmesh_ignore().

590 {
591  _coords[0] = xx;
592 
593 #if LIBMESH_DIM == 2
594  _coords[1] = xy;
595  _coords[2] = yx;
596  _coords[3] = yy;
597 #elif LIBMESH_DIM == 1
598  libmesh_assert_equal_to (xy, 0);
599  libmesh_assert_equal_to (yx, 0);
600  libmesh_assert_equal_to (yy, 0);
601  libmesh_ignore(xy, yx, yy);
602 #endif
603 
604 #if LIBMESH_DIM == 3
605  _coords[1] = xy;
606  _coords[2] = xz;
607  _coords[3] = yx;
608  _coords[4] = yy;
609  _coords[5] = yz;
610  _coords[6] = zx;
611  _coords[7] = zy;
612  _coords[8] = zz;
613 #else
614  libmesh_assert_equal_to (xz, 0);
615  libmesh_assert_equal_to (yz, 0);
616  libmesh_assert_equal_to (zx, 0);
617  libmesh_assert_equal_to (zy, 0);
618  libmesh_assert_equal_to (zz, 0);
619  libmesh_ignore(xz, yz, zx, zy, zz);
620 #endif
621 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444
void libmesh_ignore(const Args &...)

◆ TypeTensor() [4/7]

template<typename T >
template<typename T2 >
libMesh::TypeTensor< T >::TypeTensor ( const TypeVector< T2 > &  vx)
protected

Constructor.

Assigns each vector to a different row of the tensor. We're in LIBMESH_DIM space dimensions and so LIBMESH_DIM many vectors are needed.

Definition at line 638 of file type_tensor.h.

639 {
640  libmesh_assert_equal_to (LIBMESH_DIM, 1);
641  _coords[0] = vx(0);
642 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ TypeTensor() [5/7]

template<typename T >
template<typename T2 >
libMesh::TypeTensor< T >::TypeTensor ( const TypeVector< T2 > &  vx,
const TypeVector< T2 > &  vy 
)
protected

Definition at line 646 of file type_tensor.h.

References libMesh::libmesh_ignore().

648 {
649  libmesh_assert_equal_to (LIBMESH_DIM, 2);
650 #if LIBMESH_DIM > 1
651  _coords[0] = vx(0);
652  _coords[1] = vx(1);
653  _coords[2] = vy(0);
654  _coords[3] = vy(1);
655 #else
656  libmesh_ignore(vx, vy);
657 #endif
658 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444
void libmesh_ignore(const Args &...)

◆ TypeTensor() [6/7]

template<typename T >
template<typename T2 >
libMesh::TypeTensor< T >::TypeTensor ( const TypeVector< T2 > &  vx,
const TypeVector< T2 > &  vy,
const TypeVector< T2 > &  vz 
)
protected

Definition at line 662 of file type_tensor.h.

References libMesh::libmesh_ignore().

665 {
666  libmesh_assert_equal_to (LIBMESH_DIM, 3);
667 #if LIBMESH_DIM > 2
668  _coords[0] = vx(0);
669  _coords[1] = vx(1);
670  _coords[2] = vx(2);
671  _coords[3] = vy(0);
672  _coords[4] = vy(1);
673  _coords[5] = vy(2);
674  _coords[6] = vz(0);
675  _coords[7] = vz(1);
676  _coords[8] = vz(2);
677 #else
678  libmesh_ignore(vx, vy, vz);
679 #endif
680 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444
void libmesh_ignore(const Args &...)

◆ TypeTensor() [7/7]

template<typename T >
template<typename T2 >
libMesh::TypeTensor< T >::TypeTensor ( const TypeTensor< T2 > &  p)
inline

Copy-constructor.

Definition at line 628 of file type_tensor.h.

References libMesh::TypeTensor< T >::_coords.

629 {
630  // copy the nodes from vector p to me
631  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
632  _coords[i] = p._coords[i];
633 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ ~TypeTensor()

template<typename T >
libMesh::TypeTensor< T >::~TypeTensor ( )
inline

Destructor.

Definition at line 687 of file type_tensor.h.

688 {
689 }

Member Function Documentation

◆ add()

template<typename T >
template<typename T2 >
void libMesh::TypeTensor< T >::add ( const TypeTensor< T2 > &  p)
inline

Add to this tensor without creating a temporary.

Definition at line 840 of file type_tensor.h.

References libMesh::TypeTensor< T >::_coords.

841 {
842  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
843  _coords[i] += p._coords[i];
844 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ add_scaled()

template<typename T >
template<typename T2 >
void libMesh::TypeTensor< T >::add_scaled ( const TypeTensor< T2 > &  p,
const T &  factor 
)
inline

Add a scaled tensor to this tensor without creating a temporary.

Definition at line 851 of file type_tensor.h.

References libMesh::TypeTensor< T >::_coords.

Referenced by libMesh::HPCoarsenTest::add_projection(), LinearElasticityWithContact::compute_stresses(), LargeDeformationElasticity::compute_stresses(), libMesh::MeshFunction::hessian(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::System::point_hessian(), and libMesh::HPCoarsenTest::select_refinement().

852 {
853  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
854  _coords[i] += factor*p._coords[i];
855 
856 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ assign()

template<typename T >
template<typename T2 >
void libMesh::TypeTensor< T >::assign ( const TypeTensor< T2 > &  p)
inline

Assign to this tensor without creating a temporary.

Definition at line 696 of file type_tensor.h.

References libMesh::TypeTensor< T >::_coords.

697 {
698  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
699  _coords[i] = p._coords[i];
700 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ column()

template<typename T >
TypeVector< T > libMesh::TypeTensor< T >::column ( const unsigned int  r) const
inline
Returns
A copy of one column of the tensor as a TypeVector.

Definition at line 779 of file type_tensor.h.

References libMesh::TypeVector< T >::_coords.

Referenced by TypeTensorTest::testRowCol().

780 {
781  TypeVector<T> return_vector;
782 
783  for (unsigned int i=0; i<LIBMESH_DIM; i++)
784  return_vector._coords[i] = _coords[i*LIBMESH_DIM + r];
785 
786  return return_vector;
787 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ contract()

template<typename T >
template<typename T2 >
CompareTypes< T, T2 >::supertype libMesh::TypeTensor< T >::contract ( const TypeTensor< T2 > &  t) const
inline

Multiply 2 tensors together to return a scalar, i.e.

Multiply 2 tensors together, i.e.

\( \sum_{ij} A_{ij} B_{ij} \) The tensors may contain different numeric types. Also known as the "double inner product" or "double dot product" of tensors.

Returns
The scalar-valued result, this tensor is unchanged.

sum Aij*Bij. The tensors may be of different types.

Definition at line 1268 of file type_tensor.h.

References libMesh::TypeTensor< T >::_coords.

Referenced by libMesh::HPCoarsenTest::add_projection(), libMesh::VariationalSmootherSystem::element_time_derivative(), libMesh::TensorTools::inner_product(), and libMesh::HPCoarsenTest::select_refinement().

1269 {
1270  typename CompareTypes<T,T2>::supertype sum = 0.;
1271  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1272  sum += _coords[i]*t._coords[i];
1273  return sum;
1274 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ det()

template<typename T >
T libMesh::TypeTensor< T >::det ( ) const
inline
Returns
The determinant of the tensor.

Because these are 3x3 tensors at most, we don't do an LU decomposition like DenseMatrix does.

Definition at line 1298 of file type_tensor.h.

Referenced by LargeDeformationElasticity::compute_stresses(), libMesh::VariationalSmootherSystem::element_time_derivative(), libMesh::Hex::quality(), and libMesh::Sphere::Sphere().

1299 {
1300 #if LIBMESH_DIM == 1
1301  return _coords[0];
1302 #endif
1303 
1304 #if LIBMESH_DIM == 2
1305  return (_coords[0] * _coords[3]
1306  - _coords[1] * _coords[2]);
1307 #endif
1308 
1309 #if LIBMESH_DIM == 3
1310  return (_coords[0] * _coords[4] * _coords[8]
1311  + _coords[1] * _coords[5] * _coords[6]
1312  + _coords[2] * _coords[3] * _coords[7]
1313  - _coords[0] * _coords[5] * _coords[7]
1314  - _coords[1] * _coords[3] * _coords[8]
1315  - _coords[2] * _coords[4] * _coords[6]);
1316 #endif
1317 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ inverse()

template<typename T >
TypeTensor< T > libMesh::TypeTensor< T >::inverse ( ) const
inline
Returns
The inverse of this tensor as an independent object.

Definition at line 1080 of file type_tensor.h.

Referenced by NonlinearNeoHookeCurrentConfig::calculate_stress(), libMesh::VariationalSmootherSystem::element_time_derivative(), and TypeTensorTest::testInverse().

1081 {
1082 #if LIBMESH_DIM == 1
1083  if (_coords[0] == static_cast<T>(0.))
1084  libmesh_convergence_failure();
1085  return TypeTensor(1. / _coords[0]);
1086 #endif
1087 
1088 #if LIBMESH_DIM == 2
1089  // Get convenient reference to this.
1090  const TypeTensor<T> & A = *this;
1091 
1092  // Use temporary variables, avoid multiple accesses
1093  T a = A(0,0), b = A(0,1),
1094  c = A(1,0), d = A(1,1);
1095 
1096  // Make sure det = ad - bc is not zero
1097  T my_det = a*d - b*c;
1098 
1099  if (my_det == static_cast<T>(0.))
1100  libmesh_convergence_failure();
1101 
1102  return TypeTensor(d/my_det, -b/my_det, -c/my_det, a/my_det);
1103 #endif
1104 
1105 #if LIBMESH_DIM == 3
1106  // Get convenient reference to this.
1107  const TypeTensor<T> & A = *this;
1108 
1109  T a11 = A(0,0), a12 = A(0,1), a13 = A(0,2),
1110  a21 = A(1,0), a22 = A(1,1), a23 = A(1,2),
1111  a31 = A(2,0), a32 = A(2,1), a33 = A(2,2);
1112 
1113  T my_det = a11*(a33*a22-a32*a23) - a21*(a33*a12-a32*a13) + a31*(a23*a12-a22*a13);
1114 
1115  if (my_det == static_cast<T>(0.))
1116  libmesh_convergence_failure();
1117 
1118  // Inline comment characters are for lining up columns.
1119  return TypeTensor( (a33*a22-a32*a23)/my_det, -(a33*a12-a32*a13)/my_det, (a23*a12-a22*a13)/my_det,
1120  -(a33*a21-a31*a23)/my_det, (a33*a11-a31*a13)/my_det, -(a23*a11-a21*a13)/my_det,
1121  (a32*a21-a31*a22)/my_det, -(a32*a11-a31*a12)/my_det, (a22*a11-a21*a12)/my_det);
1122 #endif
1123 }
TypeTensor()
Empty constructor.
Definition: type_tensor.h:510
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ is_zero()

template<typename T >
bool libMesh::TypeTensor< T >::is_zero ( ) const
inline
Returns
True if all values in the tensor are zero

Definition at line 1288 of file type_tensor.h.

Referenced by TypeTensorTest::testIsZero().

1289 {
1290  for (const auto & val : _coords)
1291  if (val != T(0))
1292  return false;
1293  return true;
1294 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ left_multiply()

template<typename T >
template<typename T2 >
TypeVector< typename CompareTypes< T, T2 >::supertype > libMesh::TypeTensor< T >::left_multiply ( const TypeVector< T2 > &  p) const
inline

Left-multiply this tensor by a vector, i.e.

matrix-vector product. The tensor and vector may contain different numeric types.

Returns
A copy of the result vector, this tensor is unchanged.

Definition at line 1211 of file type_tensor.h.

Referenced by libMesh::operator*().

1212 {
1213  TypeVector<typename CompareTypes<T,T2>::supertype> returnval;
1214  for (unsigned int i=0; i<LIBMESH_DIM; i++)
1215  for (unsigned int j=0; j<LIBMESH_DIM; j++)
1216  returnval(i) += p(j)*(*this)(j,i);
1217 
1218  return returnval;
1219 }

◆ norm()

template<typename T >
auto libMesh::TypeTensor< T >::norm ( ) const -> decltype(std::norm(T()))
inline
Returns
The Frobenius norm of the tensor, i.e. the square-root of the sum of the elements squared.

Definition at line 1280 of file type_tensor.h.

References libMesh::TensorTools::norm_sq().

1281 {
1282  return std::sqrt(this->norm_sq());
1283 }
auto norm_sq() const -> decltype(std::norm(T()))
Definition: type_tensor.h:1348

◆ norm_sq()

template<typename T >
auto libMesh::TypeTensor< T >::norm_sq ( ) const -> decltype(std::norm(T()))
inline
Returns
The Frobenius norm of the tensor squared, i.e. sum of the element magnitudes squared.

Definition at line 1348 of file type_tensor.h.

References libMesh::TensorTools::norm_sq(), and libMesh::Real.

Referenced by libMesh::UniformRefinementEstimator::_estimate_error(), libMesh::ExactErrorEstimator::find_squared_element_error(), and libMesh::HPCoarsenTest::select_refinement().

1349 {
1350  Real sum = 0.;
1351  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1352  sum += TensorTools::norm_sq(_coords[i]);
1353  return sum;
1354 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
auto norm_sq(const T &a) -> decltype(std::norm(a))
Definition: tensor_tools.h:104

◆ operator()() [1/2]

template<typename T >
const T & libMesh::TypeTensor< T >::operator() ( const unsigned int  i,
const unsigned int  j 
) const
inline
Returns
A const reference to the (i,j) entry of the tensor.

Definition at line 706 of file type_tensor.h.

708 {
709  libmesh_assert_less (i, 3);
710  libmesh_assert_less (j, 3);
711 
712 #if LIBMESH_DIM < 3
713  const static T my_zero = 0;
714  if (i >= LIBMESH_DIM || j >= LIBMESH_DIM)
715  return my_zero;
716 #endif
717 
718  return _coords[i*LIBMESH_DIM+j];
719 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ operator()() [2/2]

template<typename T >
T & libMesh::TypeTensor< T >::operator() ( const unsigned int  i,
const unsigned int  j 
)
inline
Returns
A writable reference to the (i,j) entry of the tensor.

Definition at line 725 of file type_tensor.h.

727 {
728 #if LIBMESH_DIM < 3
729 
730  libmesh_error_msg_if(i >= LIBMESH_DIM || j >= LIBMESH_DIM,
731  "ERROR: You are assigning to a tensor component that is out of range for the compiled LIBMESH_DIM!");
732 
733 #endif
734 
735  libmesh_assert_less (i, LIBMESH_DIM);
736  libmesh_assert_less (j, LIBMESH_DIM);
737 
738  return _coords[i*LIBMESH_DIM+j];
739 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ operator*() [1/3]

template<typename T >
template<typename Scalar >
auto libMesh::TypeTensor< T >::operator* ( const Scalar &  scalar) const -> typename boostcopy::enable_if_c< ScalarTraits<Scalar>::value, TypeTensor<decltype(T() * scalar)>>::type
inline

Multiply this tensor by a scalar value.

Returns
A copy of the result, this tensor is unchanged.

Definition at line 965 of file type_tensor.h.

968 {
969  typedef decltype((*this)(0, 0) * factor) TS;
970 
971 
972 #if LIBMESH_DIM == 1
973  return TypeTensor<TS>(_coords[0]*factor);
974 #endif
975 
976 #if LIBMESH_DIM == 2
977  return TypeTensor<TS>(_coords[0]*factor,
978  _coords[1]*factor,
979  _coords[2]*factor,
980  _coords[3]*factor);
981 #endif
982 
983 #if LIBMESH_DIM == 3
984  return TypeTensor<TS>(_coords[0]*factor,
985  _coords[1]*factor,
986  _coords[2]*factor,
987  _coords[3]*factor,
988  _coords[4]*factor,
989  _coords[5]*factor,
990  _coords[6]*factor,
991  _coords[7]*factor,
992  _coords[8]*factor);
993 #endif
994 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ operator*() [2/3]

template<typename T >
template<typename T2 >
TypeTensor< typename CompareTypes< T, T2 >::supertype > libMesh::TypeTensor< T >::operator* ( const TypeTensor< T2 > &  p) const
inline

Multiply 2 tensors together, i.e.

matrix-matrix product. The tensors may contain different numeric types.

Returns
A copy of the result, this tensor is unchanged.

Definition at line 1233 of file type_tensor.h.

1234 {
1235  TypeTensor<typename CompareTypes<T, T2>::supertype> returnval;
1236  for (unsigned int i=0; i<LIBMESH_DIM; i++)
1237  for (unsigned int j=0; j<LIBMESH_DIM; j++)
1238  for (unsigned int k=0; k<LIBMESH_DIM; k++)
1239  returnval(i,j) += (*this)(i,k)*p(k,j);
1240 
1241  return returnval;
1242 }

◆ operator*() [3/3]

template<typename T >
template<typename T2 >
TypeVector< typename CompareTypes< T, T2 >::supertype > libMesh::TypeTensor< T >::operator* ( const TypeVector< T2 > &  p) const
inline

Right-multiply this tensor by a vector, i.e.

matrix-vector product. The tensor and vector may contain different numeric types.

Returns
A copy of the result vector, this tensor is unchanged.

Definition at line 1197 of file type_tensor.h.

1198 {
1199  TypeVector<typename CompareTypes<T,T2>::supertype> returnval;
1200  for (unsigned int i=0; i<LIBMESH_DIM; i++)
1201  for (unsigned int j=0; j<LIBMESH_DIM; j++)
1202  returnval(i) += (*this)(i,j)*p(j);
1203 
1204  return returnval;
1205 }

◆ operator*=() [1/2]

template<typename T>
template<typename Scalar , typename boostcopy::enable_if_c< ScalarTraits< Scalar >::value, int >::type = 0>
const TypeTensor<T>& libMesh::TypeTensor< T >::operator*= ( const Scalar &  factor)
inline

Multiply this tensor by a scalar value in place.

Returns
A reference to *this.

Definition at line 268 of file type_tensor.h.

References libMesh::TypeTensor< T >::_coords.

269  {
270  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
271  _coords[i] *= factor;
272 
273  return *this;
274  }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ operator*=() [2/2]

template<typename T >
template<typename T2 >
const TypeTensor< T > & libMesh::TypeTensor< T >::operator*= ( const TypeTensor< T2 > &  p)
inline

Multiply this tensor by a tensor value in place.

Returns
A reference to *this

Definition at line 1247 of file type_tensor.h.

1248 {
1249  TypeTensor<T> temp;
1250  for (unsigned int i=0; i<LIBMESH_DIM; i++)
1251  for (unsigned int j=0; j<LIBMESH_DIM; j++)
1252  for (unsigned int k=0; k<LIBMESH_DIM; k++)
1253  temp(i,j) += (*this)(i,k)*p(k,j);
1254 
1255  this->assign(temp);
1256  return *this;
1257 }
void assign(const TypeTensor< T2 > &)
Assign to this tensor without creating a temporary.
Definition: type_tensor.h:696

◆ operator+()

template<typename T >
template<typename T2 >
TypeTensor< typename CompareTypes< T, T2 >::supertype > libMesh::TypeTensor< T >::operator+ ( const TypeTensor< T2 > &  p) const
inline

Add another tensor to this tensor.

Returns
A copy of the result, this tensor is unchanged.

Definition at line 794 of file type_tensor.h.

References libMesh::TypeTensor< T >::_coords.

795 {
796 
797 #if LIBMESH_DIM == 1
798  return TypeTensor<typename CompareTypes<T, T2>::supertype>(_coords[0] + p._coords[0]);
799 #endif
800 
801 #if LIBMESH_DIM == 2
802  return TypeTensor<typename CompareTypes<T, T2>::supertype>(_coords[0] + p._coords[0],
803  _coords[1] + p._coords[1],
804  0.,
805  _coords[2] + p._coords[2],
806  _coords[3] + p._coords[3]);
807 #endif
808 
809 #if LIBMESH_DIM == 3
810  return TypeTensor<typename CompareTypes<T, T2>::supertype>(_coords[0] + p._coords[0],
811  _coords[1] + p._coords[1],
812  _coords[2] + p._coords[2],
813  _coords[3] + p._coords[3],
814  _coords[4] + p._coords[4],
815  _coords[5] + p._coords[5],
816  _coords[6] + p._coords[6],
817  _coords[7] + p._coords[7],
818  _coords[8] + p._coords[8]);
819 #endif
820 
821 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ operator+=()

template<typename T >
template<typename T2 >
const TypeTensor< T > & libMesh::TypeTensor< T >::operator+= ( const TypeTensor< T2 > &  p)
inline

Add to this tensor.

Returns
A reference to *this.

Definition at line 828 of file type_tensor.h.

829 {
830  this->add (p);
831 
832  return *this;
833 }
void add(const TypeTensor< T2 > &)
Add to this tensor without creating a temporary.
Definition: type_tensor.h:840

◆ operator-() [1/2]

template<typename T >
template<typename T2 >
TypeTensor< typename CompareTypes< T, T2 >::supertype > libMesh::TypeTensor< T >::operator- ( const TypeTensor< T2 > &  p) const
inline

Subtract a tensor from this tensor.

Returns
A copy of the result, this tensor is unchanged.

Definition at line 864 of file type_tensor.h.

References libMesh::TypeTensor< T >::_coords.

865 {
866 
867 #if LIBMESH_DIM == 1
868  return TypeTensor<typename CompareTypes<T, T2>::supertype>(_coords[0] - p._coords[0]);
869 #endif
870 
871 #if LIBMESH_DIM == 2
872  return TypeTensor<typename CompareTypes<T, T2>::supertype>(_coords[0] - p._coords[0],
873  _coords[1] - p._coords[1],
874  0.,
875  _coords[2] - p._coords[2],
876  _coords[3] - p._coords[3]);
877 #endif
878 
879 #if LIBMESH_DIM == 3
880  return TypeTensor<typename CompareTypes<T, T2>::supertype>(_coords[0] - p._coords[0],
881  _coords[1] - p._coords[1],
882  _coords[2] - p._coords[2],
883  _coords[3] - p._coords[3],
884  _coords[4] - p._coords[4],
885  _coords[5] - p._coords[5],
886  _coords[6] - p._coords[6],
887  _coords[7] - p._coords[7],
888  _coords[8] - p._coords[8]);
889 #endif
890 
891 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ operator-() [2/2]

template<typename T >
TypeTensor< T > libMesh::TypeTensor< T >::operator- ( ) const
inline
Returns
The negative of this tensor in a separate copy.

Definition at line 931 of file type_tensor.h.

932 {
933 
934 #if LIBMESH_DIM == 1
935  return TypeTensor(-_coords[0]);
936 #endif
937 
938 #if LIBMESH_DIM == 2
939  return TypeTensor(-_coords[0],
940  -_coords[1],
941  -_coords[2],
942  -_coords[3]);
943 #endif
944 
945 #if LIBMESH_DIM == 3
946  return TypeTensor(-_coords[0],
947  -_coords[1],
948  -_coords[2],
949  -_coords[3],
950  -_coords[4],
951  -_coords[5],
952  -_coords[6],
953  -_coords[7],
954  -_coords[8]);
955 #endif
956 
957 }
TypeTensor()
Empty constructor.
Definition: type_tensor.h:510
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ operator-=()

template<typename T >
template<typename T2 >
const TypeTensor< T > & libMesh::TypeTensor< T >::operator-= ( const TypeTensor< T2 > &  p)
inline

Subtract from this tensor.

Returns
A reference to *this.

Definition at line 898 of file type_tensor.h.

899 {
900  this->subtract (p);
901 
902  return *this;
903 }
void subtract(const TypeTensor< T2 > &)
Subtract from this tensor without creating a temporary.
Definition: type_tensor.h:910

◆ operator/()

template<typename T >
template<typename Scalar >
boostcopy::enable_if_c< ScalarTraits< Scalar >::value, TypeTensor< typename CompareTypes< T, Scalar >::supertype > >::type libMesh::TypeTensor< T >::operator/ ( const Scalar &  factor) const
inline

Divide each entry of this tensor by a scalar value.

Returns
A copy of the result, this tensor is unchanged.

Definition at line 1015 of file type_tensor.h.

1016 {
1017  libmesh_assert_not_equal_to (factor, static_cast<T>(0.));
1018 
1019  typedef typename CompareTypes<T, Scalar>::supertype TS;
1020 
1021 #if LIBMESH_DIM == 1
1022  return TypeTensor<TS>(_coords[0]/factor);
1023 #endif
1024 
1025 #if LIBMESH_DIM == 2
1026  return TypeTensor<TS>(_coords[0]/factor,
1027  _coords[1]/factor,
1028  _coords[2]/factor,
1029  _coords[3]/factor);
1030 #endif
1031 
1032 #if LIBMESH_DIM == 3
1033  return TypeTensor<TS>(_coords[0]/factor,
1034  _coords[1]/factor,
1035  _coords[2]/factor,
1036  _coords[3]/factor,
1037  _coords[4]/factor,
1038  _coords[5]/factor,
1039  _coords[6]/factor,
1040  _coords[7]/factor,
1041  _coords[8]/factor);
1042 #endif
1043 
1044 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ operator/=()

template<typename T >
const TypeTensor< T > & libMesh::TypeTensor< T >::operator/= ( const T &  factor)
inline

Divide each entry of this tensor by a scalar value.

Returns
A reference to *this.

Definition at line 1180 of file type_tensor.h.

1181 {
1182  libmesh_assert_not_equal_to (factor, static_cast<T>(0.));
1183 
1184  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1185  _coords[i] /= factor;
1186 
1187  return *this;
1188 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ operator<() [1/3]

template<>
bool libMesh::TypeTensor< Real >::operator< ( const TypeTensor< Real > &  rhs) const

Definition at line 69 of file type_tensor.C.

70 {
71  for (unsigned int i=0; i<LIBMESH_DIM; i++)
72  for (unsigned int j=0; j<LIBMESH_DIM; j++)
73  {
74  if ((*this)(i,j) < rhs(i,j))
75  return true;
76  if ((*this)(i,j) > rhs(i,j))
77  return false;
78  }
79  return false;
80 }

◆ operator<() [2/3]

template<>
bool libMesh::TypeTensor< Complex >::operator< ( const TypeTensor< Complex > &  rhs) const

Definition at line 102 of file type_tensor.C.

References std::imag(), and std::real().

103 {
104  for (unsigned int i=0; i<LIBMESH_DIM; i++)
105  for (unsigned int j=0; j<LIBMESH_DIM; j++)
106  {
107  if ((*this)(i,j).real() < rhs(i,j).real())
108  return true;
109  if ((*this)(i,j).real() > rhs(i,j).real())
110  return false;
111  if ((*this)(i,j).imag() < rhs(i,j).imag())
112  return true;
113  if ((*this)(i,j).imag() > rhs(i,j).imag())
114  return false;
115  }
116  return false;
117 }
boost::multiprecision::float128 real(const boost::multiprecision::float128 in)
boost::multiprecision::float128 imag(const boost::multiprecision::float128)

◆ operator<() [3/3]

template<typename T>
bool libMesh::TypeTensor< T >::operator< ( const TypeTensor< T > &  rhs) const
Returns
true if this tensor is "less" than another. Useful for sorting.

◆ operator=()

template<typename T>
template<typename Scalar >
boostcopy::enable_if_c< ScalarTraits<Scalar>::value, TypeTensor &>::type libMesh::TypeTensor< T >::operator= ( const Scalar &  libmesh_dbg_varp)
inline

Assignment-from-scalar operator.

Used only to zero out vectors.

Returns
A reference to *this.

Definition at line 153 of file type_tensor.h.

References libMesh::TypeTensor< T >::zero().

154  { libmesh_assert_equal_to (p, Scalar(0)); this->zero(); return *this; }
void zero()
Set all entries of the tensor to 0.
Definition: type_tensor.h:1338

◆ operator==()

template<typename T >
bool libMesh::TypeTensor< T >::operator== ( const TypeTensor< T > &  rhs) const
inline
Returns
true if two tensors are equal, false otherwise.

Definition at line 1360 of file type_tensor.h.

References libMesh::TypeTensor< T >::_coords, and libMesh::TOLERANCE.

1361 {
1362 #if LIBMESH_DIM == 1
1363  return (std::abs(_coords[0] - rhs._coords[0])
1364  < TOLERANCE);
1365 #endif
1366 
1367 #if LIBMESH_DIM == 2
1368  return ((std::abs(_coords[0] - rhs._coords[0]) +
1369  std::abs(_coords[1] - rhs._coords[1]) +
1370  std::abs(_coords[2] - rhs._coords[2]) +
1371  std::abs(_coords[3] - rhs._coords[3]))
1372  < 4.*TOLERANCE);
1373 #endif
1374 
1375 #if LIBMESH_DIM == 3
1376  return ((std::abs(_coords[0] - rhs._coords[0]) +
1377  std::abs(_coords[1] - rhs._coords[1]) +
1378  std::abs(_coords[2] - rhs._coords[2]) +
1379  std::abs(_coords[3] - rhs._coords[3]) +
1380  std::abs(_coords[4] - rhs._coords[4]) +
1381  std::abs(_coords[5] - rhs._coords[5]) +
1382  std::abs(_coords[6] - rhs._coords[6]) +
1383  std::abs(_coords[7] - rhs._coords[7]) +
1384  std::abs(_coords[8] - rhs._coords[8]))
1385  < 9.*TOLERANCE);
1386 #endif
1387 
1388 }
static constexpr Real TOLERANCE
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ operator>() [1/3]

template<>
bool libMesh::TypeTensor< Real >::operator> ( const TypeTensor< Real > &  rhs) const

Definition at line 85 of file type_tensor.C.

86 {
87  for (unsigned int i=0; i<LIBMESH_DIM; i++)
88  for (unsigned int j=0; j<LIBMESH_DIM; j++)
89  {
90  if ((*this)(i,j) > rhs(i,j))
91  return true;
92  if ((*this)(i,j) < rhs(i,j))
93  return false;
94  }
95  return false;
96 }

◆ operator>() [2/3]

template<>
bool libMesh::TypeTensor< Complex >::operator> ( const TypeTensor< Complex > &  rhs) const

Definition at line 122 of file type_tensor.C.

References std::imag(), and std::real().

123 {
124  for (unsigned int i=0; i<LIBMESH_DIM; i++)
125  for (unsigned int j=0; j<LIBMESH_DIM; j++)
126  {
127  if ((*this)(i,j).real() > rhs(i,j).real())
128  return true;
129  if ((*this)(i,j).real() < rhs(i,j).real())
130  return false;
131  if ((*this)(i,j).imag() > rhs(i,j).imag())
132  return true;
133  if ((*this)(i,j).imag() < rhs(i,j).imag())
134  return false;
135  }
136  return false;
137 }
boost::multiprecision::float128 real(const boost::multiprecision::float128 in)
boost::multiprecision::float128 imag(const boost::multiprecision::float128)

◆ operator>() [3/3]

template<typename T>
bool libMesh::TypeTensor< T >::operator> ( const TypeTensor< T > &  rhs) const
Returns
true if this tensor is "greater" than another.

◆ print()

template<typename T >
void libMesh::TypeTensor< T >::print ( std::ostream &  os = libMesh::out) const

Formatted print to a stream which defaults to libMesh::out.

Definition at line 1391 of file type_tensor.h.

1392 {
1393 #if LIBMESH_DIM == 1
1394 
1395  os << "x=" << (*this)(0,0) << std::endl;
1396 
1397 #endif
1398 #if LIBMESH_DIM == 2
1399 
1400  os << "(xx,xy)=("
1401  << std::setw(8) << (*this)(0,0) << ", "
1402  << std::setw(8) << (*this)(0,1) << ")"
1403  << std::endl;
1404  os << "(yx,yy)=("
1405  << std::setw(8) << (*this)(1,0) << ", "
1406  << std::setw(8) << (*this)(1,1) << ")"
1407  << std::endl;
1408 
1409 #endif
1410 #if LIBMESH_DIM == 3
1411 
1412  os << "(xx,xy,xz)=("
1413  << std::setw(8) << (*this)(0,0) << ", "
1414  << std::setw(8) << (*this)(0,1) << ", "
1415  << std::setw(8) << (*this)(0,2) << ")"
1416  << std::endl;
1417  os << "(yx,yy,yz)=("
1418  << std::setw(8) << (*this)(1,0) << ", "
1419  << std::setw(8) << (*this)(1,1) << ", "
1420  << std::setw(8) << (*this)(1,2) << ")"
1421  << std::endl;
1422  os << "(zx,zy,zz)=("
1423  << std::setw(8) << (*this)(2,0) << ", "
1424  << std::setw(8) << (*this)(2,1) << ", "
1425  << std::setw(8) << (*this)(2,2) << ")"
1426  << std::endl;
1427 #endif
1428 }

◆ row()

template<typename T >
TypeVector< T > libMesh::TypeTensor< T >::row ( const unsigned int  r) const
inline
Returns
A copy of one row of the tensor as a TypeVector.

Definition at line 765 of file type_tensor.h.

References libMesh::TypeVector< T >::_coords.

Referenced by TypeTensorTest::testRowCol().

766 {
767  TypeVector<T> return_vector;
768 
769  for (unsigned int j=0; j<LIBMESH_DIM; j++)
770  return_vector._coords[j] = _coords[r*LIBMESH_DIM + j];
771 
772  return return_vector;
773 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ slice() [1/2]

template<typename T >
ConstTypeTensorColumn< T > libMesh::TypeTensor< T >::slice ( const unsigned int  i) const
inline
Returns
A proxy for the \( i^{th} \) column of the tensor.

Definition at line 745 of file type_tensor.h.

746 {
747  libmesh_assert_less (i, LIBMESH_DIM);
748  return ConstTypeTensorColumn<T>(*this, i);
749 }

◆ slice() [2/2]

template<typename T >
TypeTensorColumn< T > libMesh::TypeTensor< T >::slice ( const unsigned int  i)
inline
Returns
A writable proxy for the \( i^{th} \) column of the tensor.

Definition at line 755 of file type_tensor.h.

756 {
757  libmesh_assert_less (i, LIBMESH_DIM);
758  return TypeTensorColumn<T>(*this, i);
759 }

◆ solve()

template<typename T >
void libMesh::TypeTensor< T >::solve ( const TypeVector< T > &  b,
TypeVector< T > &  x 
) const
inline

Solve the 2x2 or 3x3 system of equations A*x = b for x by directly inverting A and multiplying it by b.

Returns
The solution in the x vector.

Definition at line 1129 of file type_tensor.h.

Referenced by InfFERadialTest::base_point(), and libMesh::InfFEMap::inverse_map().

1130 {
1131 #if LIBMESH_DIM == 1
1132  if (_coords[0] == static_cast<T>(0.))
1133  libmesh_convergence_failure();
1134  x(0) = b(0) / _coords[0];
1135 #endif
1136 
1137 #if LIBMESH_DIM == 2
1138  T my_det = _coords[0]*_coords[3] - _coords[1]*_coords[2];
1139 
1140  if (my_det == static_cast<T>(0.))
1141  libmesh_convergence_failure();
1142 
1143  T my_det_inv = 1./my_det;
1144 
1145  x(0) = my_det_inv*( _coords[3]*b(0) - _coords[1]*b(1));
1146  x(1) = my_det_inv*(-_coords[2]*b(0) + _coords[0]*b(1));
1147 #endif
1148 
1149 #if LIBMESH_DIM == 3
1150  T my_det =
1151  // a11*(a33 *a22 - a32 *a23)
1152  _coords[0]*(_coords[8]*_coords[4] - _coords[7]*_coords[5])
1153  // -a21*(a33 *a12 - a32 *a13)
1154  -_coords[3]*(_coords[8]*_coords[1] - _coords[7]*_coords[2]) +
1155  // +a31*(a23 *a12 - a22 *a13)
1156  +_coords[6]*(_coords[5]*_coords[1] - _coords[4]*_coords[2]);
1157 
1158  if (my_det == static_cast<T>(0.))
1159  libmesh_convergence_failure();
1160 
1161  T my_det_inv = 1./my_det;
1162  x(0) = my_det_inv*((_coords[8]*_coords[4] - _coords[7]*_coords[5])*b(0) -
1163  (_coords[8]*_coords[1] - _coords[7]*_coords[2])*b(1) +
1164  (_coords[5]*_coords[1] - _coords[4]*_coords[2])*b(2));
1165 
1166  x(1) = my_det_inv*((_coords[6]*_coords[5] - _coords[8]*_coords[3])*b(0) +
1167  (_coords[8]*_coords[0] - _coords[6]*_coords[2])*b(1) -
1168  (_coords[5]*_coords[0] - _coords[3]*_coords[2])*b(2));
1169 
1170  x(2) = my_det_inv*((_coords[7]*_coords[3] - _coords[6]*_coords[4])*b(0) -
1171  (_coords[7]*_coords[0] - _coords[6]*_coords[1])*b(1) +
1172  (_coords[4]*_coords[0] - _coords[3]*_coords[1])*b(2));
1173 #endif
1174 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ subtract()

template<typename T >
template<typename T2 >
void libMesh::TypeTensor< T >::subtract ( const TypeTensor< T2 > &  p)
inline

Subtract from this tensor without creating a temporary.

Definition at line 910 of file type_tensor.h.

References libMesh::TypeTensor< T >::_coords.

911 {
912  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
913  _coords[i] -= p._coords[i];
914 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ subtract_scaled()

template<typename T >
template<typename T2 >
void libMesh::TypeTensor< T >::subtract_scaled ( const TypeTensor< T2 > &  p,
const T &  factor 
)
inline

Subtract a scaled value from this tensor without creating a temporary.

Definition at line 921 of file type_tensor.h.

References libMesh::TypeTensor< T >::_coords.

Referenced by libMesh::HPCoarsenTest::select_refinement().

922 {
923  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
924  _coords[i] -= factor*p._coords[i];
925 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ tr()

template<typename T >
T libMesh::TypeTensor< T >::tr ( ) const
inline
Returns
The trace of the tensor.

Definition at line 1321 of file type_tensor.h.

1322 {
1323 #if LIBMESH_DIM == 1
1324  return _coords[0];
1325 #endif
1326 
1327 #if LIBMESH_DIM == 2
1328  return _coords[0] + _coords[3];
1329 #endif
1330 
1331 #if LIBMESH_DIM == 3
1332  return _coords[0] + _coords[4] + _coords[8];
1333 #endif
1334 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ transpose()

template<typename T >
TypeTensor< T > libMesh::TypeTensor< T >::transpose ( ) const
inline
Returns
The transpose of this tensor (with complex numbers not conjugated).

Definition at line 1050 of file type_tensor.h.

Referenced by assemble_shell(), LargeDeformationElasticity::compute_stresses(), libMesh::VariationalSmootherSystem::element_time_derivative(), GradDivExactSolution::grad(), and CurlCurlExactSolution::grad().

1051 {
1052 #if LIBMESH_DIM == 1
1053  return TypeTensor(_coords[0]);
1054 #endif
1055 
1056 #if LIBMESH_DIM == 2
1057  return TypeTensor(_coords[0],
1058  _coords[2],
1059  _coords[1],
1060  _coords[3]);
1061 #endif
1062 
1063 #if LIBMESH_DIM == 3
1064  return TypeTensor(_coords[0],
1065  _coords[3],
1066  _coords[6],
1067  _coords[1],
1068  _coords[4],
1069  _coords[7],
1070  _coords[2],
1071  _coords[5],
1072  _coords[8]);
1073 #endif
1074 }
TypeTensor()
Empty constructor.
Definition: type_tensor.h:510
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

◆ write_unformatted()

template<typename T >
void libMesh::TypeTensor< T >::write_unformatted ( std::ostream &  out_stream,
const bool  newline = true 
) const

Unformatted print to the stream out.

Simply prints the elements of the tensor separated by spaces and newlines.

Definition at line 39 of file type_tensor.C.

References libMesh::libmesh_assert().

41 {
42  libmesh_assert (out_stream);
43 
44  out_stream << std::setiosflags(std::ios::showpoint)
45  << (*this)(0,0) << " "
46  << (*this)(0,1) << " "
47  << (*this)(0,2) << " ";
48  if (newline)
49  out_stream << '\n';
50 
51  out_stream << std::setiosflags(std::ios::showpoint)
52  << (*this)(1,0) << " "
53  << (*this)(1,1) << " "
54  << (*this)(1,2) << " ";
55  if (newline)
56  out_stream << '\n';
57 
58  out_stream << std::setiosflags(std::ios::showpoint)
59  << (*this)(2,0) << " "
60  << (*this)(2,1) << " "
61  << (*this)(2,2) << " ";
62  if (newline)
63  out_stream << '\n';
64 }
libmesh_assert(ctx)

◆ zero()

template<typename T >
void libMesh::TypeTensor< T >::zero ( )
inline

Set all entries of the tensor to 0.

Definition at line 1338 of file type_tensor.h.

Referenced by LinearElasticityWithContact::compute_stresses(), LargeDeformationElasticity::compute_stresses(), NonlinearNeoHookeCurrentConfig::init_for_qp(), libMesh::TensorValue< T >::operator=(), and libMesh::TypeTensor< T >::operator=().

1339 {
1340  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1341  _coords[i] = 0.;
1342 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:444

Friends And Related Function Documentation

◆ operator<<

template<typename T>
std::ostream& operator<< ( std::ostream &  os,
const TypeTensor< T > &  t 
)
friend

Formatted print as above but supports the syntax:

std::cout << t << std::endl;

Definition at line 426 of file type_tensor.h.

427  {
428  t.print(os);
429  return os;
430  }

◆ TypeTensor

template<typename T>
template<typename T2 >
friend class TypeTensor
friend

Definition at line 55 of file type_tensor.h.

Member Data Documentation

◆ _coords

template<typename T>
T libMesh::TypeTensor< T >::_coords[LIBMESH_DIM *LIBMESH_DIM]
protected

The documentation for this class was generated from the following files: