libMesh::TensorValue< T > Class Template Reference

#include <exact_solution.h>

Inheritance diagram for libMesh::TensorValue< T >:

Public Types

typedef T value_type
 
typedef std::tuple< unsigned int, unsigned int > index_type
 

Public Member Functions

 TensorValue ()
 
 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)
 
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)
 
template<typename T2 >
 TensorValue (const TypeVector< T2 > &vx)
 
template<typename T2 >
 TensorValue (const TypeVector< T2 > &vx, const TypeVector< T2 > &vy)
 
template<typename T2 >
 TensorValue (const TypeVector< T2 > &vx, const TypeVector< T2 > &vy, const TypeVector< T2 > &vz)
 
template<typename T2 >
 TensorValue (const TensorValue< T2 > &p)
 
template<typename T2 >
 TensorValue (const TypeTensor< T2 > &p)
 
 TensorValue (const TypeTensor< Real > &p_re, const TypeTensor< Real > &p_im)
 
template<typename Scalar >
boostcopy::enable_if_c< ScalarTraits< Scalar >::value, TensorValue & >::type operator= (const Scalar &libmesh_dbg_var(p))
 
template<typename T2 >
void assign (const TypeTensor< T2 > &)
 
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
 
template<typename T2 >
const TypeTensor< T > & operator+= (const TypeTensor< T2 > &)
 
template<typename T2 >
void add (const TypeTensor< T2 > &)
 
template<typename T2 >
void add_scaled (const TypeTensor< T2 > &, const T)
 
template<typename T2 >
TypeTensor< typename CompareTypes< T, T2 >::supertype > operator- (const TypeTensor< T2 > &) const
 
TypeTensor< T > operator- () const
 
template<typename T2 >
const TypeTensor< T > & operator-= (const TypeTensor< T2 > &)
 
template<typename T2 >
void subtract (const TypeTensor< T2 > &)
 
template<typename T2 >
void subtract_scaled (const TypeTensor< T2 > &, const T)
 
template<typename Scalar >
auto operator* (const Scalar scalar) const -> typename boostcopy::enable_if_c< ScalarTraits< Scalar >::value, TypeTensor< decltype(T() *scalar)>>::type
 
template<typename T2 >
TypeTensor< T > operator* (const TypeTensor< T2 > &) const
 
template<typename T2 >
TypeVector< typename CompareTypes< T, T2 >::supertype > operator* (const TypeVector< T2 > &) const
 
template<typename Scalar , typename boostcopy::enable_if_c< ScalarTraits< Scalar >::value, int >::type = 0>
const TypeTensor< T > & operator*= (const Scalar factor)
 
template<typename T2 >
const TypeTensor< T > & operator*= (const TypeTensor< T2 > &)
 
template<typename Scalar >
boostcopy::enable_if_c< ScalarTraits< Scalar >::value, TypeTensor< typename CompareTypes< T, Scalar >::supertype > >::type operator/ (const Scalar) const
 
const TypeTensor< T > & operator/= (const T)
 
template<typename T2 >
CompareTypes< T, T2 >::supertype contract (const TypeTensor< T2 > &) const
 
TypeTensor< T > transpose () const
 
TypeTensor< T > inverse () const
 
void solve (const TypeVector< T > &b, TypeVector< T > &x) const
 
Real size () const
 
Real norm () const
 
Real size_sq () const
 
Real norm_sq () const
 
det () const
 
tr () const
 
void zero ()
 
bool operator== (const TypeTensor< T > &rhs) const
 
bool operator< (const TypeTensor< T > &rhs) const
 
template<>
bool operator< (const TypeTensor< Real > &rhs) const
 
template<>
bool operator< (const TypeTensor< Complex > &rhs) const
 
bool operator> (const TypeTensor< T > &rhs) const
 
template<>
bool operator> (const TypeTensor< Real > &rhs) const
 
template<>
bool operator> (const TypeTensor< Complex > &rhs) const
 
void print (std::ostream &os=libMesh::out) const
 
void write_unformatted (std::ostream &out, const bool newline=true) const
 

Protected Attributes

_coords [LIBMESH_DIM *LIBMESH_DIM]
 

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 ( )
inline

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 
)
inlineexplicit

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 
)
inlineexplicit

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)
inline

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 
)
inline

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 
)
inline

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)
inline

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)
inline

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 
)
inline

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)
inlineinherited

Add to this tensor without creating a temporary.

Definition at line 797 of file type_tensor.h.

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

798 {
799  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
800  _coords[i] += p._coords[i];
801 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
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 
)
inlineinherited

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

Definition at line 808 of file type_tensor.h.

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

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

809 {
810  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
811  _coords[i] += factor*p._coords[i];
812 
813 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:438

◆ assign()

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

Assign to this tensor without creating a temporary.

Definition at line 668 of file type_tensor.h.

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

669 {
670  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
671  _coords[i] = p._coords[i];
672 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
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
inlineinherited

Multiply 2 tensors together to return a scalar, 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.

Multiply 2 tensors together, i.e. sum Aij*Bij. The tensors may be of different types.

Definition at line 1204 of file type_tensor.h.

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

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

1205 {
1206  typename CompareTypes<T,T2>::supertype sum = 0.;
1207  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1208  sum += _coords[i]*t._coords[i];
1209  return sum;
1210 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:438

◆ det()

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

Referenced by libMesh::Sphere::Sphere().

1238 {
1239 #if LIBMESH_DIM == 1
1240  return _coords[0];
1241 #endif
1242 
1243 #if LIBMESH_DIM == 2
1244  return (_coords[0] * _coords[3]
1245  - _coords[1] * _coords[2]);
1246 #endif
1247 
1248 #if LIBMESH_DIM == 3
1249  return (_coords[0] * _coords[4] * _coords[8]
1250  + _coords[1] * _coords[5] * _coords[6]
1251  + _coords[2] * _coords[3] * _coords[7]
1252  - _coords[0] * _coords[5] * _coords[7]
1253  - _coords[1] * _coords[3] * _coords[8]
1254  - _coords[2] * _coords[4] * _coords[6]);
1255 #endif
1256 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:438

◆ inverse()

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

Definition at line 1037 of file type_tensor.h.

References A.

1038 {
1039 #if LIBMESH_DIM == 1
1040  if (_coords[0] == static_cast<T>(0.))
1041  libmesh_convergence_failure();
1042  return TypeTensor(1. / _coords[0]);
1043 #endif
1044 
1045 #if LIBMESH_DIM == 2
1046  // Get convenient reference to this.
1047  const TypeTensor<T> & A = *this;
1048 
1049  // Use temporary variables, avoid multiple accesses
1050  T a = A(0,0), b = A(0,1),
1051  c = A(1,0), d = A(1,1);
1052 
1053  // Make sure det = ad - bc is not zero
1054  T my_det = a*d - b*c;
1055 
1056  if (my_det == static_cast<T>(0.))
1057  libmesh_convergence_failure();
1058 
1059  return TypeTensor(d/my_det, -b/my_det, -c/my_det, a/my_det);
1060 #endif
1061 
1062 #if LIBMESH_DIM == 3
1063  // Get convenient reference to this.
1064  const TypeTensor<T> & A = *this;
1065 
1066  T a11 = A(0,0), a12 = A(0,1), a13 = A(0,2),
1067  a21 = A(1,0), a22 = A(1,1), a23 = A(1,2),
1068  a31 = A(2,0), a32 = A(2,1), a33 = A(2,2);
1069 
1070  T my_det = a11*(a33*a22-a32*a23) - a21*(a33*a12-a32*a13) + a31*(a23*a12-a22*a13);
1071 
1072  if (my_det == static_cast<T>(0.))
1073  libmesh_convergence_failure();
1074 
1075  // Inline comment characters are for lining up columns.
1076  return TypeTensor( (a33*a22-a32*a23)/my_det, -(a33*a12-a32*a13)/my_det, (a23*a12-a22*a13)/my_det,
1077  -(a33*a21-a31*a23)/my_det, (a33*a11-a31*a13)/my_det, -(a23*a11-a21*a13)/my_det,
1078  (a32*a21-a31*a22)/my_det, -(a32*a11-a31*a12)/my_det, (a22*a11-a21*a12)/my_det);
1079 #endif
1080 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:438
static PetscErrorCode Mat * A

◆ norm()

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

Definition at line 1228 of file type_tensor.h.

References libMesh::TensorTools::norm_sq().

Referenced by libMesh::System::calculate_norm().

1229 {
1230  return std::sqrt(this->norm_sq());
1231 }
Real norm_sq() const
Definition: type_tensor.h:1299

◆ norm_sq()

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

Definition at line 1299 of file type_tensor.h.

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

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

1300 {
1301  Real sum = 0.;
1302  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1303  sum += TensorTools::norm_sq(_coords[i]);
1304  return sum;
1305 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
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
inlineinherited
Returns
A const reference to the (i,j) entry of the tensor.

Definition at line 678 of file type_tensor.h.

680 {
681  libmesh_assert_less (i, 3);
682  libmesh_assert_less (j, 3);
683 
684 #if LIBMESH_DIM < 3
685  const static T my_zero = 0;
686  if (i >= LIBMESH_DIM || j >= LIBMESH_DIM)
687  return my_zero;
688 #endif
689 
690  return _coords[i*LIBMESH_DIM+j];
691 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:438

◆ operator()() [2/2]

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

Definition at line 697 of file type_tensor.h.

699 {
700 #if LIBMESH_DIM < 3
701 
702  if (i >= LIBMESH_DIM || j >= LIBMESH_DIM)
703  libmesh_error_msg("ERROR: You are assigning to a tensor component that is out of range for the compiled LIBMESH_DIM!");
704 
705 #endif
706 
707  libmesh_assert_less (i, LIBMESH_DIM);
708  libmesh_assert_less (j, LIBMESH_DIM);
709 
710  return _coords[i*LIBMESH_DIM+j];
711 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
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
inlineinherited

Multiply this tensor by a scalar value.

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

Definition at line 922 of file type_tensor.h.

925 {
926  typedef decltype((*this)(0, 0) * factor) TS;
927 
928 
929 #if LIBMESH_DIM == 1
930  return TypeTensor<TS>(_coords[0]*factor);
931 #endif
932 
933 #if LIBMESH_DIM == 2
934  return TypeTensor<TS>(_coords[0]*factor,
935  _coords[1]*factor,
936  _coords[2]*factor,
937  _coords[3]*factor);
938 #endif
939 
940 #if LIBMESH_DIM == 3
941  return TypeTensor<TS>(_coords[0]*factor,
942  _coords[1]*factor,
943  _coords[2]*factor,
944  _coords[3]*factor,
945  _coords[4]*factor,
946  _coords[5]*factor,
947  _coords[6]*factor,
948  _coords[7]*factor,
949  _coords[8]*factor);
950 #endif
951 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
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
inlineinherited

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

1170 {
1171  TypeTensor<T> returnval;
1172  for (unsigned int i=0; i<LIBMESH_DIM; i++)
1173  for (unsigned int j=0; j<LIBMESH_DIM; j++)
1174  for (unsigned int k=0; k<LIBMESH_DIM; k++)
1175  returnval(i,j) += (*this)(i,k)*p(k,j);
1176 
1177  return returnval;
1178 }

◆ operator*() [3/3]

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

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

1155 {
1156  TypeVector<typename CompareTypes<T,T2>::supertype> returnval;
1157  for (unsigned int i=0; i<LIBMESH_DIM; i++)
1158  for (unsigned int j=0; j<LIBMESH_DIM; j++)
1159  returnval(i) += (*this)(i,j)*p(j);
1160 
1161  return returnval;
1162 }

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

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]
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)
inlineinherited

Multiply this tensor by a tensor value in place

Returns
A reference to *this

Definition at line 1183 of file type_tensor.h.

1184 {
1185  TypeTensor<T> temp;
1186  for (unsigned int i=0; i<LIBMESH_DIM; i++)
1187  for (unsigned int j=0; j<LIBMESH_DIM; j++)
1188  for (unsigned int k=0; k<LIBMESH_DIM; k++)
1189  temp(i,j) += (*this)(i,k)*p(k,j);
1190 
1191  this->assign(temp);
1192  return *this;
1193 }
void assign(const TypeTensor< T2 > &)
Definition: type_tensor.h:668

◆ operator+()

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

Add another tensor to this tensor.

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

Definition at line 751 of file type_tensor.h.

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

752 {
753 
754 #if LIBMESH_DIM == 1
755  return TypeTensor(_coords[0] + p._coords[0]);
756 #endif
757 
758 #if LIBMESH_DIM == 2
759  return TypeTensor(_coords[0] + p._coords[0],
760  _coords[1] + p._coords[1],
761  0.,
762  _coords[2] + p._coords[2],
763  _coords[3] + p._coords[3]);
764 #endif
765 
766 #if LIBMESH_DIM == 3
767  return TypeTensor(_coords[0] + p._coords[0],
768  _coords[1] + p._coords[1],
769  _coords[2] + p._coords[2],
770  _coords[3] + p._coords[3],
771  _coords[4] + p._coords[4],
772  _coords[5] + p._coords[5],
773  _coords[6] + p._coords[6],
774  _coords[7] + p._coords[7],
775  _coords[8] + p._coords[8]);
776 #endif
777 
778 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:438

◆ operator+=()

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

Add to this tensor.

Returns
A reference to *this.

Definition at line 785 of file type_tensor.h.

786 {
787  this->add (p);
788 
789  return *this;
790 }
void add(const TypeTensor< T2 > &)
Definition: type_tensor.h:797

◆ operator-() [1/2]

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

Subtract a tensor from this tensor.

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

Definition at line 821 of file type_tensor.h.

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

822 {
823 
824 #if LIBMESH_DIM == 1
825  return TypeTensor(_coords[0] - p._coords[0]);
826 #endif
827 
828 #if LIBMESH_DIM == 2
829  return TypeTensor(_coords[0] - p._coords[0],
830  _coords[1] - p._coords[1],
831  0.,
832  _coords[2] - p._coords[2],
833  _coords[3] - p._coords[3]);
834 #endif
835 
836 #if LIBMESH_DIM == 3
837  return TypeTensor(_coords[0] - p._coords[0],
838  _coords[1] - p._coords[1],
839  _coords[2] - p._coords[2],
840  _coords[3] - p._coords[3],
841  _coords[4] - p._coords[4],
842  _coords[5] - p._coords[5],
843  _coords[6] - p._coords[6],
844  _coords[7] - p._coords[7],
845  _coords[8] - p._coords[8]);
846 #endif
847 
848 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:438

◆ operator-() [2/2]

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

Definition at line 888 of file type_tensor.h.

889 {
890 
891 #if LIBMESH_DIM == 1
892  return TypeTensor(-_coords[0]);
893 #endif
894 
895 #if LIBMESH_DIM == 2
896  return TypeTensor(-_coords[0],
897  -_coords[1],
898  -_coords[2],
899  -_coords[3]);
900 #endif
901 
902 #if LIBMESH_DIM == 3
903  return TypeTensor(-_coords[0],
904  -_coords[1],
905  -_coords[2],
906  -_coords[3],
907  -_coords[4],
908  -_coords[5],
909  -_coords[6],
910  -_coords[7],
911  -_coords[8]);
912 #endif
913 
914 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:438

◆ operator-=()

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

Subtract from this tensor.

Returns
A reference to *this.

Definition at line 855 of file type_tensor.h.

856 {
857  this->subtract (p);
858 
859  return *this;
860 }
void subtract(const TypeTensor< T2 > &)
Definition: type_tensor.h:867

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

Divide each entry of this tensor by a scalar value.

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

Definition at line 972 of file type_tensor.h.

973 {
974  libmesh_assert_not_equal_to (factor, static_cast<T>(0.));
975 
976  typedef typename CompareTypes<T, Scalar>::supertype TS;
977 
978 #if LIBMESH_DIM == 1
979  return TypeTensor<TS>(_coords[0]/factor);
980 #endif
981 
982 #if LIBMESH_DIM == 2
983  return TypeTensor<TS>(_coords[0]/factor,
984  _coords[1]/factor,
985  _coords[2]/factor,
986  _coords[3]/factor);
987 #endif
988 
989 #if LIBMESH_DIM == 3
990  return TypeTensor<TS>(_coords[0]/factor,
991  _coords[1]/factor,
992  _coords[2]/factor,
993  _coords[3]/factor,
994  _coords[4]/factor,
995  _coords[5]/factor,
996  _coords[6]/factor,
997  _coords[7]/factor,
998  _coords[8]/factor);
999 #endif
1000 
1001 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:438

◆ operator/=()

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

Divide each entry of this tensor by a scalar value.

Returns
A reference to *this.

Definition at line 1137 of file type_tensor.h.

1138 {
1139  libmesh_assert_not_equal_to (factor, static_cast<T>(0.));
1140 
1141  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1142  _coords[i] /= factor;
1143 
1144  return *this;
1145 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:438

◆ operator<() [1/3]

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

Definition at line 113 of file type_tensor.C.

114 {
115  for (unsigned int i=0; i<LIBMESH_DIM; i++)
116  for (unsigned int j=0; j<LIBMESH_DIM; j++)
117  {
118  if ((*this)(i,j) < rhs(i,j))
119  return true;
120  if ((*this)(i,j) > rhs(i,j))
121  return false;
122  }
123  return false;
124 }

◆ operator<() [2/3]

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

Definition at line 146 of file type_tensor.C.

147 {
148  for (unsigned int i=0; i<LIBMESH_DIM; i++)
149  for (unsigned int j=0; j<LIBMESH_DIM; j++)
150  {
151  if ((*this)(i,j).real() < rhs(i,j).real())
152  return true;
153  if ((*this)(i,j).real() > rhs(i,j).real())
154  return false;
155  if ((*this)(i,j).imag() < rhs(i,j).imag())
156  return true;
157  if ((*this)(i,j).imag() > rhs(i,j).imag())
158  return false;
159  }
160  return false;
161 }

◆ operator<() [3/3]

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)
inline

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; }

◆ operator==()

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

Definition at line 1311 of file type_tensor.h.

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

1312 {
1313 #if LIBMESH_DIM == 1
1314  return (std::abs(_coords[0] - rhs._coords[0])
1315  < TOLERANCE);
1316 #endif
1317 
1318 #if LIBMESH_DIM == 2
1319  return ((std::abs(_coords[0] - rhs._coords[0]) +
1320  std::abs(_coords[1] - rhs._coords[1]) +
1321  std::abs(_coords[2] - rhs._coords[2]) +
1322  std::abs(_coords[3] - rhs._coords[3]))
1323  < 4.*TOLERANCE);
1324 #endif
1325 
1326 #if LIBMESH_DIM == 3
1327  return ((std::abs(_coords[0] - rhs._coords[0]) +
1328  std::abs(_coords[1] - rhs._coords[1]) +
1329  std::abs(_coords[2] - rhs._coords[2]) +
1330  std::abs(_coords[3] - rhs._coords[3]) +
1331  std::abs(_coords[4] - rhs._coords[4]) +
1332  std::abs(_coords[5] - rhs._coords[5]) +
1333  std::abs(_coords[6] - rhs._coords[6]) +
1334  std::abs(_coords[7] - rhs._coords[7]) +
1335  std::abs(_coords[8] - rhs._coords[8]))
1336  < 9.*TOLERANCE);
1337 #endif
1338 
1339 }
double abs(double a)
static const Real TOLERANCE
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:438

◆ operator>() [1/3]

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

Definition at line 129 of file type_tensor.C.

130 {
131  for (unsigned int i=0; i<LIBMESH_DIM; i++)
132  for (unsigned int j=0; j<LIBMESH_DIM; j++)
133  {
134  if ((*this)(i,j) > rhs(i,j))
135  return true;
136  if ((*this)(i,j) < rhs(i,j))
137  return false;
138  }
139  return false;
140 }

◆ operator>() [2/3]

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

Definition at line 166 of file type_tensor.C.

167 {
168  for (unsigned int i=0; i<LIBMESH_DIM; i++)
169  for (unsigned int j=0; j<LIBMESH_DIM; j++)
170  {
171  if ((*this)(i,j).real() > rhs(i,j).real())
172  return true;
173  if ((*this)(i,j).real() < rhs(i,j).real())
174  return false;
175  if ((*this)(i,j).imag() > rhs(i,j).imag())
176  return true;
177  if ((*this)(i,j).imag() < rhs(i,j).imag())
178  return false;
179  }
180  return false;
181 }

◆ operator>() [3/3]

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.

Definition at line 39 of file type_tensor.C.

40 {
41 #if LIBMESH_DIM == 1
42 
43  os << "x=" << (*this)(0,0) << std::endl;
44 
45 #endif
46 #if LIBMESH_DIM == 2
47 
48  os << "(xx,xy)=("
49  << std::setw(8) << (*this)(0,0) << ", "
50  << std::setw(8) << (*this)(0,1) << ")"
51  << std::endl;
52  os << "(yx,yy)=("
53  << std::setw(8) << (*this)(1,0) << ", "
54  << std::setw(8) << (*this)(1,1) << ")"
55  << std::endl;
56 
57 #endif
58 #if LIBMESH_DIM == 3
59 
60  os << "(xx,xy,xz)=("
61  << std::setw(8) << (*this)(0,0) << ", "
62  << std::setw(8) << (*this)(0,1) << ", "
63  << std::setw(8) << (*this)(0,2) << ")"
64  << std::endl;
65  os << "(yx,yy,yz)=("
66  << std::setw(8) << (*this)(1,0) << ", "
67  << std::setw(8) << (*this)(1,1) << ", "
68  << std::setw(8) << (*this)(1,2) << ")"
69  << std::endl;
70  os << "(zx,zy,zz)=("
71  << std::setw(8) << (*this)(2,0) << ", "
72  << std::setw(8) << (*this)(2,1) << ", "
73  << std::setw(8) << (*this)(2,2) << ")"
74  << std::endl;
75 #endif
76 }

◆ row()

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

Definition at line 737 of file type_tensor.h.

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

738 {
739  TypeVector<T> return_vector;
740 
741  for (unsigned int j=0; j<LIBMESH_DIM; j++)
742  return_vector._coords[j] = _coords[r*LIBMESH_DIM + j];
743 
744  return return_vector;
745 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:438

◆ size()

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

1218 {
1219  libmesh_deprecated();
1220  return this->norm();
1221 }
Real norm() const
Definition: type_tensor.h:1228

◆ size_sq()

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

References libMesh::TensorTools::norm_sq().

1289 {
1290  libmesh_deprecated();
1291  return this->norm_sq();
1292 }
Real norm_sq() const
Definition: type_tensor.h:1299

◆ slice() [1/2]

template<typename T >
ConstTypeTensorColumn< T > libMesh::TypeTensor< T >::slice ( const unsigned int  i) const
inlineinherited
Returns
A 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 ConstTypeTensorColumn<T>(*this, i);
721 }

◆ slice() [2/2]

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

Definition at line 727 of file type_tensor.h.

728 {
729  libmesh_assert_less (i, LIBMESH_DIM);
730  return TypeTensorColumn<T>(*this, i);
731 }

◆ solve()

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

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

1087 {
1088 #if LIBMESH_DIM == 1
1089  if (_coords[0] == static_cast<T>(0.))
1090  libmesh_convergence_failure();
1091  x(0) = b(0) / _coords[0];
1092 #endif
1093 
1094 #if LIBMESH_DIM == 2
1095  T my_det = _coords[0]*_coords[3] - _coords[1]*_coords[2];
1096 
1097  if (my_det == static_cast<T>(0.))
1098  libmesh_convergence_failure();
1099 
1100  T my_det_inv = 1./my_det;
1101 
1102  x(0) = my_det_inv*( _coords[3]*b(0) - _coords[1]*b(1));
1103  x(1) = my_det_inv*(-_coords[2]*b(0) + _coords[0]*b(1));
1104 #endif
1105 
1106 #if LIBMESH_DIM == 3
1107  T my_det =
1108  // a11*(a33 *a22 - a32 *a23)
1109  _coords[0]*(_coords[8]*_coords[4] - _coords[7]*_coords[5])
1110  // -a21*(a33 *a12 - a32 *a13)
1111  -_coords[3]*(_coords[8]*_coords[1] - _coords[7]*_coords[2]) +
1112  // +a31*(a23 *a12 - a22 *a13)
1113  +_coords[6]*(_coords[5]*_coords[1] - _coords[4]*_coords[2]);
1114 
1115  if (my_det == static_cast<T>(0.))
1116  libmesh_convergence_failure();
1117 
1118  T my_det_inv = 1./my_det;
1119  x(0) = my_det_inv*((_coords[8]*_coords[4] - _coords[7]*_coords[5])*b(0) -
1120  (_coords[8]*_coords[1] - _coords[7]*_coords[2])*b(1) +
1121  (_coords[5]*_coords[1] - _coords[4]*_coords[2])*b(2));
1122 
1123  x(1) = my_det_inv*((_coords[6]*_coords[5] - _coords[8]*_coords[3])*b(0) +
1124  (_coords[8]*_coords[0] - _coords[6]*_coords[2])*b(1) -
1125  (_coords[5]*_coords[0] - _coords[3]*_coords[2])*b(2));
1126 
1127  x(2) = my_det_inv*((_coords[7]*_coords[3] - _coords[6]*_coords[4])*b(0) -
1128  (_coords[7]*_coords[0] - _coords[6]*_coords[1])*b(1) +
1129  (_coords[4]*_coords[0] - _coords[3]*_coords[1])*b(2));
1130 #endif
1131 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:438

◆ subtract()

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

Subtract from this tensor without creating a temporary.

Definition at line 867 of file type_tensor.h.

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

868 {
869  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
870  _coords[i] -= p._coords[i];
871 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
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 
)
inlineinherited

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

Definition at line 878 of file type_tensor.h.

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

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

879 {
880  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
881  _coords[i] -= factor*p._coords[i];
882 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:438

◆ tr()

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

Definition at line 1260 of file type_tensor.h.

1261 {
1262 #if LIBMESH_DIM == 1
1263  return _coords[0];
1264 #endif
1265 
1266 #if LIBMESH_DIM == 2
1267  return _coords[0] + _coords[3];
1268 #endif
1269 
1270 #if LIBMESH_DIM == 3
1271  return _coords[0] + _coords[4] + _coords[8];
1272 #endif
1273 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:438

◆ transpose()

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

Definition at line 1007 of file type_tensor.h.

1008 {
1009 #if LIBMESH_DIM == 1
1010  return TypeTensor(_coords[0]);
1011 #endif
1012 
1013 #if LIBMESH_DIM == 2
1014  return TypeTensor(_coords[0],
1015  _coords[2],
1016  _coords[1],
1017  _coords[3]);
1018 #endif
1019 
1020 #if LIBMESH_DIM == 3
1021  return TypeTensor(_coords[0],
1022  _coords[3],
1023  _coords[6],
1024  _coords[1],
1025  _coords[4],
1026  _coords[7],
1027  _coords[2],
1028  _coords[5],
1029  _coords[8]);
1030 #endif
1031 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
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.

Definition at line 83 of file type_tensor.C.

85 {
86  libmesh_assert (out_stream);
87 
88  out_stream << std::setiosflags(std::ios::showpoint)
89  << (*this)(0,0) << " "
90  << (*this)(0,1) << " "
91  << (*this)(0,2) << " ";
92  if (newline)
93  out_stream << '\n';
94 
95  out_stream << std::setiosflags(std::ios::showpoint)
96  << (*this)(1,0) << " "
97  << (*this)(1,1) << " "
98  << (*this)(1,2) << " ";
99  if (newline)
100  out_stream << '\n';
101 
102  out_stream << std::setiosflags(std::ios::showpoint)
103  << (*this)(2,0) << " "
104  << (*this)(2,1) << " "
105  << (*this)(2,2) << " ";
106  if (newline)
107  out_stream << '\n';
108 }

◆ zero()

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

Set all entries of the tensor to 0.

Definition at line 1277 of file type_tensor.h.

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

1278 {
1279  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1280  _coords[i] = 0.;
1281 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:438

Member Data Documentation

◆ _coords


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