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< T > 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...
 
Real size () const
 
Real norm () const
 
Real size_sq () const
 
Real norm_sq () 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, 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 121 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 116 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 787 of file type_tensor.h.

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

788 {
789  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
790  _coords[i] += p._coords[i];
791 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

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

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

799 {
800  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
801  _coords[i] += factor*p._coords[i];
802 
803 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

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

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

659 {
660  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
661  _coords[i] = p._coords[i];
662 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

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

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

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

1195 {
1196  typename CompareTypes<T,T2>::supertype sum = 0.;
1197  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1198  sum += _coords[i]*t._coords[i];
1199  return sum;
1200 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

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

1228 {
1229 #if LIBMESH_DIM == 1
1230  return _coords[0];
1231 #endif
1232 
1233 #if LIBMESH_DIM == 2
1234  return (_coords[0] * _coords[3]
1235  - _coords[1] * _coords[2]);
1236 #endif
1237 
1238 #if LIBMESH_DIM == 3
1239  return (_coords[0] * _coords[4] * _coords[8]
1240  + _coords[1] * _coords[5] * _coords[6]
1241  + _coords[2] * _coords[3] * _coords[7]
1242  - _coords[0] * _coords[5] * _coords[7]
1243  - _coords[1] * _coords[3] * _coords[8]
1244  - _coords[2] * _coords[4] * _coords[6]);
1245 #endif
1246 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

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

1028 {
1029 #if LIBMESH_DIM == 1
1030  if (_coords[0] == static_cast<T>(0.))
1031  libmesh_convergence_failure();
1032  return TypeTensor(1. / _coords[0]);
1033 #endif
1034 
1035 #if LIBMESH_DIM == 2
1036  // Get convenient reference to this.
1037  const TypeTensor<T> & A = *this;
1038 
1039  // Use temporary variables, avoid multiple accesses
1040  T a = A(0,0), b = A(0,1),
1041  c = A(1,0), d = A(1,1);
1042 
1043  // Make sure det = ad - bc is not zero
1044  T my_det = a*d - b*c;
1045 
1046  if (my_det == static_cast<T>(0.))
1047  libmesh_convergence_failure();
1048 
1049  return TypeTensor(d/my_det, -b/my_det, -c/my_det, a/my_det);
1050 #endif
1051 
1052 #if LIBMESH_DIM == 3
1053  // Get convenient reference to this.
1054  const TypeTensor<T> & A = *this;
1055 
1056  T a11 = A(0,0), a12 = A(0,1), a13 = A(0,2),
1057  a21 = A(1,0), a22 = A(1,1), a23 = A(1,2),
1058  a31 = A(2,0), a32 = A(2,1), a33 = A(2,2);
1059 
1060  T my_det = a11*(a33*a22-a32*a23) - a21*(a33*a12-a32*a13) + a31*(a23*a12-a22*a13);
1061 
1062  if (my_det == static_cast<T>(0.))
1063  libmesh_convergence_failure();
1064 
1065  // Inline comment characters are for lining up columns.
1066  return TypeTensor( (a33*a22-a32*a23)/my_det, -(a33*a12-a32*a13)/my_det, (a23*a12-a22*a13)/my_det,
1067  -(a33*a21-a31*a23)/my_det, (a33*a11-a31*a13)/my_det, -(a23*a11-a21*a13)/my_det,
1068  (a32*a21-a31*a22)/my_det, -(a32*a11-a31*a12)/my_det, (a22*a11-a21*a12)/my_det);
1069 #endif
1070 }
TypeTensor()
Empty constructor.
Definition: type_tensor.h:504
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

◆ norm()

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

Definition at line 1218 of file type_tensor.h.

References libMesh::TensorTools::norm_sq().

1219 {
1220  return std::sqrt(this->norm_sq());
1221 }
Real norm_sq() const
Definition: type_tensor.h:1289

◆ norm_sq()

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

Definition at line 1289 of file type_tensor.h.

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

1290 {
1291  Real sum = 0.;
1292  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1293  sum += TensorTools::norm_sq(_coords[i]);
1294  return sum;
1295 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438
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 668 of file type_tensor.h.

670 {
671  libmesh_assert_less (i, 3);
672  libmesh_assert_less (j, 3);
673 
674 #if LIBMESH_DIM < 3
675  const static T my_zero = 0;
676  if (i >= LIBMESH_DIM || j >= LIBMESH_DIM)
677  return my_zero;
678 #endif
679 
680  return _coords[i*LIBMESH_DIM+j];
681 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

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

689 {
690 #if LIBMESH_DIM < 3
691 
692  if (i >= LIBMESH_DIM || j >= LIBMESH_DIM)
693  libmesh_error_msg("ERROR: You are assigning to a tensor component that is out of range for the compiled LIBMESH_DIM!");
694 
695 #endif
696 
697  libmesh_assert_less (i, LIBMESH_DIM);
698  libmesh_assert_less (j, LIBMESH_DIM);
699 
700  return _coords[i*LIBMESH_DIM+j];
701 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

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

915 {
916  typedef decltype((*this)(0, 0) * factor) TS;
917 
918 
919 #if LIBMESH_DIM == 1
920  return TypeTensor<TS>(_coords[0]*factor);
921 #endif
922 
923 #if LIBMESH_DIM == 2
924  return TypeTensor<TS>(_coords[0]*factor,
925  _coords[1]*factor,
926  _coords[2]*factor,
927  _coords[3]*factor);
928 #endif
929 
930 #if LIBMESH_DIM == 3
931  return TypeTensor<TS>(_coords[0]*factor,
932  _coords[1]*factor,
933  _coords[2]*factor,
934  _coords[3]*factor,
935  _coords[4]*factor,
936  _coords[5]*factor,
937  _coords[6]*factor,
938  _coords[7]*factor,
939  _coords[8]*factor);
940 #endif
941 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

◆ operator*() [2/3]

template<typename T >
template<typename T2 >
TypeTensor< T > 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 1159 of file type_tensor.h.

1160 {
1161  TypeTensor<T> returnval;
1162  for (unsigned int i=0; i<LIBMESH_DIM; i++)
1163  for (unsigned int j=0; j<LIBMESH_DIM; j++)
1164  for (unsigned int k=0; k<LIBMESH_DIM; k++)
1165  returnval(i,j) += (*this)(i,k)*p(k,j);
1166 
1167  return returnval;
1168 }

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

1145 {
1146  TypeVector<typename CompareTypes<T,T2>::supertype> returnval;
1147  for (unsigned int i=0; i<LIBMESH_DIM; i++)
1148  for (unsigned int j=0; j<LIBMESH_DIM; j++)
1149  returnval(i) += (*this)(i,j)*p(j);
1150 
1151  return returnval;
1152 }

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

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

260  {
261  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
262  _coords[i] *= factor;
263 
264  return *this;
265  }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

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

1174 {
1175  TypeTensor<T> temp;
1176  for (unsigned int i=0; i<LIBMESH_DIM; i++)
1177  for (unsigned int j=0; j<LIBMESH_DIM; j++)
1178  for (unsigned int k=0; k<LIBMESH_DIM; k++)
1179  temp(i,j) += (*this)(i,k)*p(k,j);
1180 
1181  this->assign(temp);
1182  return *this;
1183 }
void assign(const TypeTensor< T2 > &)
Assign to this tensor without creating a temporary.
Definition: type_tensor.h:658

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

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

742 {
743 
744 #if LIBMESH_DIM == 1
745  return TypeTensor(_coords[0] + p._coords[0]);
746 #endif
747 
748 #if LIBMESH_DIM == 2
749  return TypeTensor(_coords[0] + p._coords[0],
750  _coords[1] + p._coords[1],
751  0.,
752  _coords[2] + p._coords[2],
753  _coords[3] + p._coords[3]);
754 #endif
755 
756 #if LIBMESH_DIM == 3
757  return TypeTensor(_coords[0] + p._coords[0],
758  _coords[1] + p._coords[1],
759  _coords[2] + p._coords[2],
760  _coords[3] + p._coords[3],
761  _coords[4] + p._coords[4],
762  _coords[5] + p._coords[5],
763  _coords[6] + p._coords[6],
764  _coords[7] + p._coords[7],
765  _coords[8] + p._coords[8]);
766 #endif
767 
768 }
TypeTensor()
Empty constructor.
Definition: type_tensor.h:504
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

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

776 {
777  this->add (p);
778 
779  return *this;
780 }
void add(const TypeTensor< T2 > &)
Add to this tensor without creating a temporary.
Definition: type_tensor.h:787

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

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

812 {
813 
814 #if LIBMESH_DIM == 1
815  return TypeTensor(_coords[0] - p._coords[0]);
816 #endif
817 
818 #if LIBMESH_DIM == 2
819  return TypeTensor(_coords[0] - p._coords[0],
820  _coords[1] - p._coords[1],
821  0.,
822  _coords[2] - p._coords[2],
823  _coords[3] - p._coords[3]);
824 #endif
825 
826 #if LIBMESH_DIM == 3
827  return TypeTensor(_coords[0] - p._coords[0],
828  _coords[1] - p._coords[1],
829  _coords[2] - p._coords[2],
830  _coords[3] - p._coords[3],
831  _coords[4] - p._coords[4],
832  _coords[5] - p._coords[5],
833  _coords[6] - p._coords[6],
834  _coords[7] - p._coords[7],
835  _coords[8] - p._coords[8]);
836 #endif
837 
838 }
TypeTensor()
Empty constructor.
Definition: type_tensor.h:504
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

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

879 {
880 
881 #if LIBMESH_DIM == 1
882  return TypeTensor(-_coords[0]);
883 #endif
884 
885 #if LIBMESH_DIM == 2
886  return TypeTensor(-_coords[0],
887  -_coords[1],
888  -_coords[2],
889  -_coords[3]);
890 #endif
891 
892 #if LIBMESH_DIM == 3
893  return TypeTensor(-_coords[0],
894  -_coords[1],
895  -_coords[2],
896  -_coords[3],
897  -_coords[4],
898  -_coords[5],
899  -_coords[6],
900  -_coords[7],
901  -_coords[8]);
902 #endif
903 
904 }
TypeTensor()
Empty constructor.
Definition: type_tensor.h:504
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

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

846 {
847  this->subtract (p);
848 
849  return *this;
850 }
void subtract(const TypeTensor< T2 > &)
Subtract from this tensor without creating a temporary.
Definition: type_tensor.h:857

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

963 {
964  libmesh_assert_not_equal_to (factor, static_cast<T>(0.));
965 
966  typedef typename CompareTypes<T, Scalar>::supertype TS;
967 
968 #if LIBMESH_DIM == 1
969  return TypeTensor<TS>(_coords[0]/factor);
970 #endif
971 
972 #if LIBMESH_DIM == 2
973  return TypeTensor<TS>(_coords[0]/factor,
974  _coords[1]/factor,
975  _coords[2]/factor,
976  _coords[3]/factor);
977 #endif
978 
979 #if LIBMESH_DIM == 3
980  return TypeTensor<TS>(_coords[0]/factor,
981  _coords[1]/factor,
982  _coords[2]/factor,
983  _coords[3]/factor,
984  _coords[4]/factor,
985  _coords[5]/factor,
986  _coords[6]/factor,
987  _coords[7]/factor,
988  _coords[8]/factor);
989 #endif
990 
991 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

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

1128 {
1129  libmesh_assert_not_equal_to (factor, static_cast<T>(0.));
1130 
1131  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1132  _coords[i] /= factor;
1133 
1134  return *this;
1135 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

◆ 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:1267

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

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

1302 {
1303 #if LIBMESH_DIM == 1
1304  return (std::abs(_coords[0] - rhs._coords[0])
1305  < TOLERANCE);
1306 #endif
1307 
1308 #if LIBMESH_DIM == 2
1309  return ((std::abs(_coords[0] - rhs._coords[0]) +
1310  std::abs(_coords[1] - rhs._coords[1]) +
1311  std::abs(_coords[2] - rhs._coords[2]) +
1312  std::abs(_coords[3] - rhs._coords[3]))
1313  < 4.*TOLERANCE);
1314 #endif
1315 
1316 #if LIBMESH_DIM == 3
1317  return ((std::abs(_coords[0] - rhs._coords[0]) +
1318  std::abs(_coords[1] - rhs._coords[1]) +
1319  std::abs(_coords[2] - rhs._coords[2]) +
1320  std::abs(_coords[3] - rhs._coords[3]) +
1321  std::abs(_coords[4] - rhs._coords[4]) +
1322  std::abs(_coords[5] - rhs._coords[5]) +
1323  std::abs(_coords[6] - rhs._coords[6]) +
1324  std::abs(_coords[7] - rhs._coords[7]) +
1325  std::abs(_coords[8] - rhs._coords[8]))
1326  < 9.*TOLERANCE);
1327 #endif
1328 
1329 }
double abs(double a)
static const Real TOLERANCE
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

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

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

728 {
729  TypeVector<T> return_vector;
730 
731  for (unsigned int j=0; j<LIBMESH_DIM; j++)
732  return_vector._coords[j] = _coords[r*LIBMESH_DIM + j];
733 
734  return return_vector;
735 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

◆ size()

template<typename T >
Real libMesh::TypeTensor< T >::size ( ) const
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 1207 of file type_tensor.h.

1208 {
1209  libmesh_deprecated();
1210  return this->norm();
1211 }
Real norm() const
Definition: type_tensor.h:1218

◆ size_sq()

template<typename T >
Real libMesh::TypeTensor< T >::size_sq ( ) const
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 1278 of file type_tensor.h.

References libMesh::TensorTools::norm_sq().

1279 {
1280  libmesh_deprecated();
1281  return this->norm_sq();
1282 }
Real norm_sq() const
Definition: type_tensor.h:1289

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

708 {
709  libmesh_assert_less (i, LIBMESH_DIM);
710  return ConstTypeTensorColumn<T>(*this, i);
711 }

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

718 {
719  libmesh_assert_less (i, LIBMESH_DIM);
720  return TypeTensorColumn<T>(*this, i);
721 }

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

1077 {
1078 #if LIBMESH_DIM == 1
1079  if (_coords[0] == static_cast<T>(0.))
1080  libmesh_convergence_failure();
1081  x(0) = b(0) / _coords[0];
1082 #endif
1083 
1084 #if LIBMESH_DIM == 2
1085  T my_det = _coords[0]*_coords[3] - _coords[1]*_coords[2];
1086 
1087  if (my_det == static_cast<T>(0.))
1088  libmesh_convergence_failure();
1089 
1090  T my_det_inv = 1./my_det;
1091 
1092  x(0) = my_det_inv*( _coords[3]*b(0) - _coords[1]*b(1));
1093  x(1) = my_det_inv*(-_coords[2]*b(0) + _coords[0]*b(1));
1094 #endif
1095 
1096 #if LIBMESH_DIM == 3
1097  T my_det =
1098  // a11*(a33 *a22 - a32 *a23)
1099  _coords[0]*(_coords[8]*_coords[4] - _coords[7]*_coords[5])
1100  // -a21*(a33 *a12 - a32 *a13)
1101  -_coords[3]*(_coords[8]*_coords[1] - _coords[7]*_coords[2]) +
1102  // +a31*(a23 *a12 - a22 *a13)
1103  +_coords[6]*(_coords[5]*_coords[1] - _coords[4]*_coords[2]);
1104 
1105  if (my_det == static_cast<T>(0.))
1106  libmesh_convergence_failure();
1107 
1108  T my_det_inv = 1./my_det;
1109  x(0) = my_det_inv*((_coords[8]*_coords[4] - _coords[7]*_coords[5])*b(0) -
1110  (_coords[8]*_coords[1] - _coords[7]*_coords[2])*b(1) +
1111  (_coords[5]*_coords[1] - _coords[4]*_coords[2])*b(2));
1112 
1113  x(1) = my_det_inv*((_coords[6]*_coords[5] - _coords[8]*_coords[3])*b(0) +
1114  (_coords[8]*_coords[0] - _coords[6]*_coords[2])*b(1) -
1115  (_coords[5]*_coords[0] - _coords[3]*_coords[2])*b(2));
1116 
1117  x(2) = my_det_inv*((_coords[7]*_coords[3] - _coords[6]*_coords[4])*b(0) -
1118  (_coords[7]*_coords[0] - _coords[6]*_coords[1])*b(1) +
1119  (_coords[4]*_coords[0] - _coords[3]*_coords[1])*b(2));
1120 #endif
1121 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

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

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

858 {
859  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
860  _coords[i] -= p._coords[i];
861 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

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

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

869 {
870  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
871  _coords[i] -= factor*p._coords[i];
872 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

◆ tr()

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

Definition at line 1250 of file type_tensor.h.

1251 {
1252 #if LIBMESH_DIM == 1
1253  return _coords[0];
1254 #endif
1255 
1256 #if LIBMESH_DIM == 2
1257  return _coords[0] + _coords[3];
1258 #endif
1259 
1260 #if LIBMESH_DIM == 3
1261  return _coords[0] + _coords[4] + _coords[8];
1262 #endif
1263 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

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

998 {
999 #if LIBMESH_DIM == 1
1000  return TypeTensor(_coords[0]);
1001 #endif
1002 
1003 #if LIBMESH_DIM == 2
1004  return TypeTensor(_coords[0],
1005  _coords[2],
1006  _coords[1],
1007  _coords[3]);
1008 #endif
1009 
1010 #if LIBMESH_DIM == 3
1011  return TypeTensor(_coords[0],
1012  _coords[3],
1013  _coords[6],
1014  _coords[1],
1015  _coords[4],
1016  _coords[7],
1017  _coords[2],
1018  _coords[5],
1019  _coords[8]);
1020 #endif
1021 }
TypeTensor()
Empty constructor.
Definition: type_tensor.h:504
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

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

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

1268 {
1269  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1270  _coords[i] = 0.;
1271 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:438

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: