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

This class defines a tensor in LIBMESH_DIM dimensional Real or Complex space. More...

#include <exact_solution.h>

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

Public Types

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

Public Member Functions

 TensorValue ()
 Empty constructor. More...
 
 TensorValue (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 >
 TensorValue (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-scalars. More...
 
template<typename T2 >
 TensorValue (const TypeVector< T2 > &vx)
 Constructor. More...
 
template<typename T2 >
 TensorValue (const TypeVector< T2 > &vx, const TypeVector< T2 > &vy)
 Constructor. More...
 
template<typename T2 >
 TensorValue (const TypeVector< T2 > &vx, const TypeVector< T2 > &vy, const TypeVector< T2 > &vz)
 Constructor. More...
 
template<typename T2 >
 TensorValue (const TensorValue< T2 > &p)
 Copy-constructor. More...
 
template<typename T2 >
 TensorValue (const TypeTensor< T2 > &p)
 Copy-constructor. More...
 
 TensorValue (const TypeTensor< Real > &p_re, const TypeTensor< Real > &p_im)
 Constructor that takes two TypeTensor<Real> representing the real and imaginary part as arguments. More...
 
template<typename Scalar >
boostcopy::enable_if_c< ScalarTraits< Scalar >::value, TensorValue & >::type operator= (const Scalar &libmesh_dbg_var(p))
 Assignment-from-scalar operator. More...
 
template<typename T2 >
void assign (const TypeTensor< T2 > &)
 Assign to this tensor without creating a temporary. 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
 
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...
 
TypeTensor< T > operator- () const
 
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...
 
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 T2 >
TypeTensor< typename CompareTypes< T, T2 >::supertype > operator* (const TypeTensor< T2 > &) const
 Multiply 2 tensors together, i.e. More...
 
template<typename T2 >
TypeVector< typename CompareTypes< T, T2 >::supertype > operator* (const TypeVector< T2 > &) const
 Multiply this tensor by a vector, i.e. 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 T2 >
const TypeTensor< T > & operator*= (const TypeTensor< T2 > &)
 Multiply this tensor by a tensor 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 >
CompareTypes< T, T2 >::supertype contract (const TypeTensor< T2 > &) const
 Multiply 2 tensors together to return a scalar, 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 size () const -> decltype(std::norm(T()))
 
auto norm () const -> decltype(std::norm(T()))
 
auto size_sq () const -> decltype(std::norm(T()))
 
auto norm_sq () const -> decltype(std::norm(T()))
 
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, const bool newline=true) const
 Unformatted print to the stream out. More...
 

Protected Attributes

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

Detailed Description

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

This class defines a tensor in LIBMESH_DIM dimensional Real or Complex space.

The typedef RealTensorValue always defines a real-valued tensor, and NumberTensorValue defines a real or complex-valued tensor depending on how the library was configured.

Author
Roy H. Stogner
Date
2004

Definition at line 52 of file exact_solution.h.

Member Typedef Documentation

◆ index_type

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

Helper typedef for generic index programming.

Definition at line 144 of file type_tensor.h.

◆ value_type

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

Helper typedef for C++98 generic programming.

Definition at line 139 of file type_tensor.h.

Constructor & Destructor Documentation

◆ TensorValue() [1/9]

template<typename T >
libMesh::TensorValue< T >::TensorValue ( )

Empty constructor.

Gives the tensor 0 in LIBMESH_DIM dimensional T space.

Definition at line 156 of file tensor_value.h.

156  :
157  TypeTensor<T> ()
158 {
159 }

◆ TensorValue() [2/9]

template<typename T >
libMesh::TensorValue< T >::TensorValue ( 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 
)
explicit

Constructor-from-T.

By default sets higher dimensional entries to 0.

Definition at line 165 of file tensor_value.h.

173  :
174  TypeTensor<T> (xx,xy,xz,yx,yy,yz,zx,zy,zz)
175 {
176 }

◆ TensorValue() [3/9]

template<typename T >
template<typename Scalar >
libMesh::TensorValue< T >::TensorValue ( 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 
)
explicit

Constructor-from-scalars.

By default sets higher dimensional entries to 0.

Definition at line 182 of file tensor_value.h.

192  :
193  TypeTensor<T> (xx,xy,xz,yx,yy,yz,zx,zy,zz)
194 {
195 }

◆ TensorValue() [4/9]

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

Constructor.

Takes 1 row vector for LIBMESH_DIM=1

Definition at line 212 of file tensor_value.h.

212  :
213  TypeTensor<T> (vx)
214 {
215 }

◆ TensorValue() [5/9]

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

Constructor.

Takes 2 row vectors for LIBMESH_DIM=2

Definition at line 222 of file tensor_value.h.

223  :
224  TypeTensor<T> (vx, vy)
225 {
226 }

◆ TensorValue() [6/9]

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

Constructor.

Takes 3 row vectors for LIBMESH_DIM=3

Definition at line 233 of file tensor_value.h.

235  :
236  TypeTensor<T> (vx, vy, vz)
237 {
238 }

◆ TensorValue() [7/9]

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

Copy-constructor.

Definition at line 202 of file tensor_value.h.

202  :
203  TypeTensor<T> (p)
204 {
205 }

◆ TensorValue() [8/9]

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

Copy-constructor.

Definition at line 245 of file tensor_value.h.

245  :
246  TypeTensor<T> (p)
247 {
248 }

◆ TensorValue() [9/9]

template<typename T >
libMesh::TensorValue< T >::TensorValue ( const TypeTensor< Real > &  p_re,
const TypeTensor< Real > &  p_im 
)

Constructor that takes two TypeTensor<Real> representing the real and imaginary part as arguments.

Definition at line 254 of file tensor_value.h.

255  :
256  TypeTensor<T> (Complex (p_re(0,0), p_im(0,0)),
257  Complex (p_re(0,1), p_im(0,1)),
258  Complex (p_re(0,2), p_im(0,2)),
259  Complex (p_re(1,0), p_im(1,0)),
260  Complex (p_re(1,1), p_im(1,1)),
261  Complex (p_re(1,2), p_im(1,2)),
262  Complex (p_re(2,0), p_im(2,0)),
263  Complex (p_re(2,1), p_im(2,1)),
264  Complex (p_re(2,2), p_im(2,2)))
265 {
266 }
std::complex< Real > Complex

Member Function Documentation

◆ add()

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

Add to this tensor without creating a temporary.

Definition at line 821 of file type_tensor.h.

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

822 {
823  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
824  _coords[i] += p._coords[i];
825 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462

◆ add_scaled()

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

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

Definition at line 832 of file type_tensor.h.

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

833 {
834  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
835  _coords[i] += factor*p._coords[i];
836 
837 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462

◆ assign()

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

Assign to this tensor without creating a temporary.

Definition at line 692 of file type_tensor.h.

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

693 {
694  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
695  _coords[i] = p._coords[i];
696 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462

◆ contract()

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

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 1229 of file type_tensor.h.

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

Referenced by libMesh::TensorTools::inner_product().

1230 {
1231  typename CompareTypes<T,T2>::supertype sum = 0.;
1232  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1233  sum += _coords[i]*t._coords[i];
1234  return sum;
1235 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462

◆ det()

template<typename T >
T libMesh::TypeTensor< T >::det ( ) const
inherited
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 1262 of file type_tensor.h.

1263 {
1264 #if LIBMESH_DIM == 1
1265  return _coords[0];
1266 #endif
1267 
1268 #if LIBMESH_DIM == 2
1269  return (_coords[0] * _coords[3]
1270  - _coords[1] * _coords[2]);
1271 #endif
1272 
1273 #if LIBMESH_DIM == 3
1274  return (_coords[0] * _coords[4] * _coords[8]
1275  + _coords[1] * _coords[5] * _coords[6]
1276  + _coords[2] * _coords[3] * _coords[7]
1277  - _coords[0] * _coords[5] * _coords[7]
1278  - _coords[1] * _coords[3] * _coords[8]
1279  - _coords[2] * _coords[4] * _coords[6]);
1280 #endif
1281 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462

◆ inverse()

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

Definition at line 1061 of file type_tensor.h.

1062 {
1063 #if LIBMESH_DIM == 1
1064  if (_coords[0] == static_cast<T>(0.))
1065  libmesh_convergence_failure();
1066  return TypeTensor(1. / _coords[0]);
1067 #endif
1068 
1069 #if LIBMESH_DIM == 2
1070  // Get convenient reference to this.
1071  const TypeTensor<T> & A = *this;
1072 
1073  // Use temporary variables, avoid multiple accesses
1074  T a = A(0,0), b = A(0,1),
1075  c = A(1,0), d = A(1,1);
1076 
1077  // Make sure det = ad - bc is not zero
1078  T my_det = a*d - b*c;
1079 
1080  if (my_det == static_cast<T>(0.))
1081  libmesh_convergence_failure();
1082 
1083  return TypeTensor(d/my_det, -b/my_det, -c/my_det, a/my_det);
1084 #endif
1085 
1086 #if LIBMESH_DIM == 3
1087  // Get convenient reference to this.
1088  const TypeTensor<T> & A = *this;
1089 
1090  T a11 = A(0,0), a12 = A(0,1), a13 = A(0,2),
1091  a21 = A(1,0), a22 = A(1,1), a23 = A(1,2),
1092  a31 = A(2,0), a32 = A(2,1), a33 = A(2,2);
1093 
1094  T my_det = a11*(a33*a22-a32*a23) - a21*(a33*a12-a32*a13) + a31*(a23*a12-a22*a13);
1095 
1096  if (my_det == static_cast<T>(0.))
1097  libmesh_convergence_failure();
1098 
1099  // Inline comment characters are for lining up columns.
1100  return TypeTensor( (a33*a22-a32*a23)/my_det, -(a33*a12-a32*a13)/my_det, (a23*a12-a22*a13)/my_det,
1101  -(a33*a21-a31*a23)/my_det, (a33*a11-a31*a13)/my_det, -(a23*a11-a21*a13)/my_det,
1102  (a32*a21-a31*a22)/my_det, -(a32*a11-a31*a12)/my_det, (a22*a11-a21*a12)/my_det);
1103 #endif
1104 }
TypeTensor()
Empty constructor.
Definition: type_tensor.h:528
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462

◆ norm()

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

Definition at line 1253 of file type_tensor.h.

References libMesh::TensorTools::norm_sq(), and std::sqrt().

1254 {
1255  return std::sqrt(this->norm_sq());
1256 }
MetaPhysicL::DualNumber< T, D > sqrt(const MetaPhysicL::DualNumber< T, D > &in)
auto norm_sq() const -> decltype(std::norm(T()))
Definition: type_tensor.h:1324

◆ norm_sq()

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

Definition at line 1324 of file type_tensor.h.

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

1325 {
1326  Real sum = 0.;
1327  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1328  sum += TensorTools::norm_sq(_coords[i]);
1329  return sum;
1330 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462
T norm_sq(std::complex< T > a)
Definition: tensor_tools.h:85
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real

◆ operator()() [1/2]

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

Definition at line 702 of file type_tensor.h.

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

◆ operator()() [2/2]

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

Definition at line 721 of file type_tensor.h.

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

◆ 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
inherited

Multiply this tensor by a scalar value.

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

Definition at line 946 of file type_tensor.h.

949 {
950  typedef decltype((*this)(0, 0) * factor) TS;
951 
952 
953 #if LIBMESH_DIM == 1
954  return TypeTensor<TS>(_coords[0]*factor);
955 #endif
956 
957 #if LIBMESH_DIM == 2
958  return TypeTensor<TS>(_coords[0]*factor,
959  _coords[1]*factor,
960  _coords[2]*factor,
961  _coords[3]*factor);
962 #endif
963 
964 #if LIBMESH_DIM == 3
965  return TypeTensor<TS>(_coords[0]*factor,
966  _coords[1]*factor,
967  _coords[2]*factor,
968  _coords[3]*factor,
969  _coords[4]*factor,
970  _coords[5]*factor,
971  _coords[6]*factor,
972  _coords[7]*factor,
973  _coords[8]*factor);
974 #endif
975 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462

◆ operator*() [2/3]

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

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 1194 of file type_tensor.h.

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

◆ operator*() [3/3]

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

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 1178 of file type_tensor.h.

1179 {
1180  TypeVector<typename CompareTypes<T,T2>::supertype> returnval;
1181  for (unsigned int i=0; i<LIBMESH_DIM; i++)
1182  for (unsigned int j=0; j<LIBMESH_DIM; j++)
1183  returnval(i) += (*this)(i,j)*p(j);
1184 
1185  return returnval;
1186 }

◆ 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)
inherited

Multiply this tensor by a scalar value in place.

Returns
A reference to *this.

Definition at line 282 of file type_tensor.h.

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

283  {
284  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
285  _coords[i] *= factor;
286 
287  return *this;
288  }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462

◆ operator*=() [2/2]

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

Multiply this tensor by a tensor value in place.

Returns
A reference to *this

Definition at line 1208 of file type_tensor.h.

1209 {
1210  TypeTensor<T> temp;
1211  for (unsigned int i=0; i<LIBMESH_DIM; i++)
1212  for (unsigned int j=0; j<LIBMESH_DIM; j++)
1213  for (unsigned int k=0; k<LIBMESH_DIM; k++)
1214  temp(i,j) += (*this)(i,k)*p(k,j);
1215 
1216  this->assign(temp);
1217  return *this;
1218 }
void assign(const TypeTensor< T2 > &)
Assign to this tensor without creating a temporary.
Definition: type_tensor.h:692

◆ operator+()

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

Add another tensor to this tensor.

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

Definition at line 775 of file type_tensor.h.

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

776 {
777 
778 #if LIBMESH_DIM == 1
779  return TypeTensor<typename CompareTypes<T, T2>::supertype>(_coords[0] + p._coords[0]);
780 #endif
781 
782 #if LIBMESH_DIM == 2
783  return TypeTensor<typename CompareTypes<T, T2>::supertype>(_coords[0] + p._coords[0],
784  _coords[1] + p._coords[1],
785  0.,
786  _coords[2] + p._coords[2],
787  _coords[3] + p._coords[3]);
788 #endif
789 
790 #if LIBMESH_DIM == 3
791  return TypeTensor<typename CompareTypes<T, T2>::supertype>(_coords[0] + p._coords[0],
792  _coords[1] + p._coords[1],
793  _coords[2] + p._coords[2],
794  _coords[3] + p._coords[3],
795  _coords[4] + p._coords[4],
796  _coords[5] + p._coords[5],
797  _coords[6] + p._coords[6],
798  _coords[7] + p._coords[7],
799  _coords[8] + p._coords[8]);
800 #endif
801 
802 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462

◆ operator+=()

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

Add to this tensor.

Returns
A reference to *this.

Definition at line 809 of file type_tensor.h.

810 {
811  this->add (p);
812 
813  return *this;
814 }
void add(const TypeTensor< T2 > &)
Add to this tensor without creating a temporary.
Definition: type_tensor.h:821

◆ operator-() [1/2]

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

Subtract a tensor from this tensor.

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

Definition at line 845 of file type_tensor.h.

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

846 {
847 
848 #if LIBMESH_DIM == 1
849  return TypeTensor<typename CompareTypes<T, T2>::supertype>(_coords[0] - p._coords[0]);
850 #endif
851 
852 #if LIBMESH_DIM == 2
853  return TypeTensor<typename CompareTypes<T, T2>::supertype>(_coords[0] - p._coords[0],
854  _coords[1] - p._coords[1],
855  0.,
856  _coords[2] - p._coords[2],
857  _coords[3] - p._coords[3]);
858 #endif
859 
860 #if LIBMESH_DIM == 3
861  return TypeTensor<typename CompareTypes<T, T2>::supertype>(_coords[0] - p._coords[0],
862  _coords[1] - p._coords[1],
863  _coords[2] - p._coords[2],
864  _coords[3] - p._coords[3],
865  _coords[4] - p._coords[4],
866  _coords[5] - p._coords[5],
867  _coords[6] - p._coords[6],
868  _coords[7] - p._coords[7],
869  _coords[8] - p._coords[8]);
870 #endif
871 
872 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462

◆ operator-() [2/2]

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

Definition at line 912 of file type_tensor.h.

913 {
914 
915 #if LIBMESH_DIM == 1
916  return TypeTensor(-_coords[0]);
917 #endif
918 
919 #if LIBMESH_DIM == 2
920  return TypeTensor(-_coords[0],
921  -_coords[1],
922  -_coords[2],
923  -_coords[3]);
924 #endif
925 
926 #if LIBMESH_DIM == 3
927  return TypeTensor(-_coords[0],
928  -_coords[1],
929  -_coords[2],
930  -_coords[3],
931  -_coords[4],
932  -_coords[5],
933  -_coords[6],
934  -_coords[7],
935  -_coords[8]);
936 #endif
937 
938 }
TypeTensor()
Empty constructor.
Definition: type_tensor.h:528
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462

◆ operator-=()

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

Subtract from this tensor.

Returns
A reference to *this.

Definition at line 879 of file type_tensor.h.

880 {
881  this->subtract (p);
882 
883  return *this;
884 }
void subtract(const TypeTensor< T2 > &)
Subtract from this tensor without creating a temporary.
Definition: type_tensor.h:891

◆ 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
inherited

Divide each entry of this tensor by a scalar value.

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

Definition at line 996 of file type_tensor.h.

997 {
998  libmesh_assert_not_equal_to (factor, static_cast<T>(0.));
999 
1000  typedef typename CompareTypes<T, Scalar>::supertype TS;
1001 
1002 #if LIBMESH_DIM == 1
1003  return TypeTensor<TS>(_coords[0]/factor);
1004 #endif
1005 
1006 #if LIBMESH_DIM == 2
1007  return TypeTensor<TS>(_coords[0]/factor,
1008  _coords[1]/factor,
1009  _coords[2]/factor,
1010  _coords[3]/factor);
1011 #endif
1012 
1013 #if LIBMESH_DIM == 3
1014  return TypeTensor<TS>(_coords[0]/factor,
1015  _coords[1]/factor,
1016  _coords[2]/factor,
1017  _coords[3]/factor,
1018  _coords[4]/factor,
1019  _coords[5]/factor,
1020  _coords[6]/factor,
1021  _coords[7]/factor,
1022  _coords[8]/factor);
1023 #endif
1024 
1025 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462

◆ operator/=()

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

Divide each entry of this tensor by a scalar value.

Returns
A reference to *this.

Definition at line 1161 of file type_tensor.h.

1162 {
1163  libmesh_assert_not_equal_to (factor, static_cast<T>(0.));
1164 
1165  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1166  _coords[i] /= factor;
1167 
1168  return *this;
1169 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462

◆ operator<()

template<typename T>
bool libMesh::TypeTensor< T >::operator< ( const TypeTensor< T > &  rhs) const
inherited
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, TensorValue &>::type libMesh::TensorValue< T >::operator= ( const Scalar &  libmesh_dbg_varp)

Assignment-from-scalar operator.

Used only to zero out tensors.

Definition at line 135 of file tensor_value.h.

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

136  { 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:1302

◆ operator==()

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

Definition at line 1336 of file type_tensor.h.

References libMesh::TypeTensor< T >::_coords, std::abs(), and libMesh::TOLERANCE.

1337 {
1338 #if LIBMESH_DIM == 1
1339  return (std::abs(_coords[0] - rhs._coords[0])
1340  < TOLERANCE);
1341 #endif
1342 
1343 #if LIBMESH_DIM == 2
1344  return ((std::abs(_coords[0] - rhs._coords[0]) +
1345  std::abs(_coords[1] - rhs._coords[1]) +
1346  std::abs(_coords[2] - rhs._coords[2]) +
1347  std::abs(_coords[3] - rhs._coords[3]))
1348  < 4.*TOLERANCE);
1349 #endif
1350 
1351 #if LIBMESH_DIM == 3
1352  return ((std::abs(_coords[0] - rhs._coords[0]) +
1353  std::abs(_coords[1] - rhs._coords[1]) +
1354  std::abs(_coords[2] - rhs._coords[2]) +
1355  std::abs(_coords[3] - rhs._coords[3]) +
1356  std::abs(_coords[4] - rhs._coords[4]) +
1357  std::abs(_coords[5] - rhs._coords[5]) +
1358  std::abs(_coords[6] - rhs._coords[6]) +
1359  std::abs(_coords[7] - rhs._coords[7]) +
1360  std::abs(_coords[8] - rhs._coords[8]))
1361  < 9.*TOLERANCE);
1362 #endif
1363 
1364 }
double abs(double a)
static const Real TOLERANCE
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462

◆ operator>()

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

◆ print()

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

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

◆ row()

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

Definition at line 761 of file type_tensor.h.

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

762 {
763  TypeVector<T> return_vector;
764 
765  for (unsigned int j=0; j<LIBMESH_DIM; j++)
766  return_vector._coords[j] = _coords[r*LIBMESH_DIM + j];
767 
768  return return_vector;
769 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462

◆ size()

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

Definition at line 1242 of file type_tensor.h.

References std::norm().

1243 {
1244  libmesh_deprecated();
1245  return this->norm();
1246 }
auto norm() const -> decltype(std::norm(T()))
Definition: type_tensor.h:1253

◆ size_sq()

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

Definition at line 1313 of file type_tensor.h.

References libMesh::TensorTools::norm_sq().

1314 {
1315  libmesh_deprecated();
1316  return this->norm_sq();
1317 }
auto norm_sq() const -> decltype(std::norm(T()))
Definition: type_tensor.h:1324

◆ slice() [1/2]

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

Definition at line 741 of file type_tensor.h.

742 {
743  libmesh_assert_less (i, LIBMESH_DIM);
744  return ConstTypeTensorColumn<T>(*this, i);
745 }

◆ slice() [2/2]

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

Definition at line 751 of file type_tensor.h.

752 {
753  libmesh_assert_less (i, LIBMESH_DIM);
754  return TypeTensorColumn<T>(*this, i);
755 }

◆ solve()

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

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 1110 of file type_tensor.h.

1111 {
1112 #if LIBMESH_DIM == 1
1113  if (_coords[0] == static_cast<T>(0.))
1114  libmesh_convergence_failure();
1115  x(0) = b(0) / _coords[0];
1116 #endif
1117 
1118 #if LIBMESH_DIM == 2
1119  T my_det = _coords[0]*_coords[3] - _coords[1]*_coords[2];
1120 
1121  if (my_det == static_cast<T>(0.))
1122  libmesh_convergence_failure();
1123 
1124  T my_det_inv = 1./my_det;
1125 
1126  x(0) = my_det_inv*( _coords[3]*b(0) - _coords[1]*b(1));
1127  x(1) = my_det_inv*(-_coords[2]*b(0) + _coords[0]*b(1));
1128 #endif
1129 
1130 #if LIBMESH_DIM == 3
1131  T my_det =
1132  // a11*(a33 *a22 - a32 *a23)
1133  _coords[0]*(_coords[8]*_coords[4] - _coords[7]*_coords[5])
1134  // -a21*(a33 *a12 - a32 *a13)
1135  -_coords[3]*(_coords[8]*_coords[1] - _coords[7]*_coords[2]) +
1136  // +a31*(a23 *a12 - a22 *a13)
1137  +_coords[6]*(_coords[5]*_coords[1] - _coords[4]*_coords[2]);
1138 
1139  if (my_det == static_cast<T>(0.))
1140  libmesh_convergence_failure();
1141 
1142  T my_det_inv = 1./my_det;
1143  x(0) = my_det_inv*((_coords[8]*_coords[4] - _coords[7]*_coords[5])*b(0) -
1144  (_coords[8]*_coords[1] - _coords[7]*_coords[2])*b(1) +
1145  (_coords[5]*_coords[1] - _coords[4]*_coords[2])*b(2));
1146 
1147  x(1) = my_det_inv*((_coords[6]*_coords[5] - _coords[8]*_coords[3])*b(0) +
1148  (_coords[8]*_coords[0] - _coords[6]*_coords[2])*b(1) -
1149  (_coords[5]*_coords[0] - _coords[3]*_coords[2])*b(2));
1150 
1151  x(2) = my_det_inv*((_coords[7]*_coords[3] - _coords[6]*_coords[4])*b(0) -
1152  (_coords[7]*_coords[0] - _coords[6]*_coords[1])*b(1) +
1153  (_coords[4]*_coords[0] - _coords[3]*_coords[1])*b(2));
1154 #endif
1155 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462

◆ subtract()

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

Subtract from this tensor without creating a temporary.

Definition at line 891 of file type_tensor.h.

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

892 {
893  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
894  _coords[i] -= p._coords[i];
895 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462

◆ subtract_scaled()

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

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

Definition at line 902 of file type_tensor.h.

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

903 {
904  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
905  _coords[i] -= factor*p._coords[i];
906 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462

◆ tr()

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

Definition at line 1285 of file type_tensor.h.

1286 {
1287 #if LIBMESH_DIM == 1
1288  return _coords[0];
1289 #endif
1290 
1291 #if LIBMESH_DIM == 2
1292  return _coords[0] + _coords[3];
1293 #endif
1294 
1295 #if LIBMESH_DIM == 3
1296  return _coords[0] + _coords[4] + _coords[8];
1297 #endif
1298 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462

◆ transpose()

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

Definition at line 1031 of file type_tensor.h.

1032 {
1033 #if LIBMESH_DIM == 1
1034  return TypeTensor(_coords[0]);
1035 #endif
1036 
1037 #if LIBMESH_DIM == 2
1038  return TypeTensor(_coords[0],
1039  _coords[2],
1040  _coords[1],
1041  _coords[3]);
1042 #endif
1043 
1044 #if LIBMESH_DIM == 3
1045  return TypeTensor(_coords[0],
1046  _coords[3],
1047  _coords[6],
1048  _coords[1],
1049  _coords[4],
1050  _coords[7],
1051  _coords[2],
1052  _coords[5],
1053  _coords[8]);
1054 #endif
1055 }
TypeTensor()
Empty constructor.
Definition: type_tensor.h:528
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462

◆ write_unformatted()

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

Unformatted print to the stream out.

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

◆ zero()

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

Set all entries of the tensor to 0.

Definition at line 1302 of file type_tensor.h.

Referenced by libMesh::TensorValue< T >::operator=(), and libMesh::TypeTensor< T >::operator=().

1303 {
1304  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1305  _coords[i] = 0.;
1306 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:462

Member Data Documentation

◆ _coords

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

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