libMesh::TypeTensor< T > Class Template Reference

#include <tensor_tools.h>

Inheritance diagram for libMesh::TypeTensor< T >:

Public Types

typedef T value_type
 

Public Member Functions

template<typename T2 >
 TypeTensor (const TypeTensor< T2 > &p)
 
 ~TypeTensor ()
 
template<typename T2 >
void assign (const TypeTensor< T2 > &)
 
template<typename Scalar >
boostcopy::enable_if_c< ScalarTraits< Scalar >::value, TypeTensor & >::type operator= (const Scalar &libmesh_dbg_var(p))
 
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
 
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)
 
TypeTensor< T > operator- () const
 
template<typename Scalar >
boostcopy::enable_if_c< ScalarTraits< Scalar >::value, TypeTensor< typename CompareTypes< T, Scalar >::supertype > >::type operator* (const Scalar) const
 
template<typename Scalar >
const TypeTensor< T > & operator*= (const Scalar)
 
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 >
TypeTensor< T > operator* (const TypeTensor< T2 > &) const
 
template<typename T2 >
CompareTypes< T, T2 >::supertype contract (const TypeTensor< T2 > &) const
 
template<typename T2 >
TypeVector< typename CompareTypes< T, T2 >::supertype > operator* (const TypeVector< 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
 
bool operator> (const TypeTensor< T > &rhs) const
 
void print (std::ostream &os=libMesh::out) const
 
void write_unformatted (std::ostream &out, const bool newline=true) const
 
template<>
bool operator< (const TypeTensor< Real > &rhs) const
 
template<>
bool operator> (const TypeTensor< Real > &rhs) const
 
template<>
bool operator< (const TypeTensor< Complex > &rhs) const
 
template<>
bool operator> (const TypeTensor< Complex > &rhs) const
 

Protected Member Functions

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

Protected Attributes

_coords [LIBMESH_DIM *LIBMESH_DIM]
 

Friends

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

Detailed Description

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

This class defines a tensor in LIBMESH_DIM dimensional space of type T. T may either be Real or Complex. The default constructor for this class is protected, suggesting that you should not instantiate one of these directly.

Author
Roy Stogner
Date
2004

Definition at line 32 of file tensor_tools.h.

Member Typedef Documentation

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

Helper typedef for C++98 generic programming.

Definition at line 117 of file type_tensor.h.

Constructor & Destructor Documentation

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

Empty constructor. Gives the tensor 0 in LIBMESH_DIM dimensions.

Definition at line 479 of file type_tensor.h.

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

Referenced by libMesh::TypeTensor< T >::inverse(), libMesh::TypeTensor< T >::operator+(), libMesh::TypeTensor< T >::operator-(), and libMesh::TypeTensor< T >::transpose().

480 {
481  _coords[0] = 0;
482 
483 #if LIBMESH_DIM > 1
484  _coords[1] = 0;
485  _coords[2] = 0;
486  _coords[3] = 0;
487 #endif
488 
489 #if LIBMESH_DIM > 2
490  _coords[4] = 0;
491  _coords[5] = 0;
492  _coords[6] = 0;
493  _coords[7] = 0;
494  _coords[8] = 0;
495 #endif
496 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
libMesh::TypeTensor< T >::TypeTensor ( const T  xx,
const T  xy = 0,
const T  xz = 0,
const T  yx = 0,
const T  yy = 0,
const T  yz = 0,
const T  zx = 0,
const T  zy = 0,
const T  zz = 0 
)
inlineexplicitprotected

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

Definition at line 502 of file type_tensor.h.

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

511 {
512  _coords[0] = xx;
513 
514 #if LIBMESH_DIM == 2
515  _coords[1] = xy;
516  _coords[2] = yx;
517  _coords[3] = yy;
518 #endif
519 
520 #if LIBMESH_DIM == 3
521  _coords[1] = xy;
522  _coords[2] = xz;
523  _coords[3] = yx;
524  _coords[4] = yy;
525  _coords[5] = yz;
526  _coords[6] = zx;
527  _coords[7] = zy;
528  _coords[8] = zz;
529 #endif
530 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
template<typename Scalar >
libMesh::TypeTensor< T >::TypeTensor ( const Scalar  xx,
const Scalar  xy = 0,
const Scalar  xz = 0,
const Scalar  yx = 0,
const Scalar  yy = 0,
const Scalar  yz = 0,
const Scalar  zx = 0,
const Scalar  zy = 0,
typename boostcopy::enable_if_c< ScalarTraits< Scalar >::value, const Scalar >::type  zz = 0 
)
inlineexplicitprotected

Constructor-from-Scalar.

Definition at line 536 of file type_tensor.h.

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

547 {
548  _coords[0] = xx;
549 
550 #if LIBMESH_DIM == 2
551  _coords[1] = xy;
552  _coords[2] = yx;
553  _coords[3] = yy;
554 #endif
555 
556 #if LIBMESH_DIM == 3
557  _coords[1] = xy;
558  _coords[2] = xz;
559  _coords[3] = yx;
560  _coords[4] = yy;
561  _coords[5] = yz;
562  _coords[6] = zx;
563  _coords[7] = zy;
564  _coords[8] = zz;
565 #endif
566 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
template<typename T2 >
libMesh::TypeTensor< T >::TypeTensor ( const TypeVector< T2 > &  vx)
protected

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

Definition at line 583 of file type_tensor.h.

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

584 {
585  libmesh_assert_equal_to (LIBMESH_DIM, 1);
586  _coords[0] = vx(0);
587 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
template<typename T2 >
libMesh::TypeTensor< T >::TypeTensor ( const TypeVector< T2 > &  vx,
const TypeVector< T2 > &  vy 
)
protected

Definition at line 591 of file type_tensor.h.

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

593 {
594  libmesh_assert_equal_to (LIBMESH_DIM, 2);
595  _coords[0] = vx(0);
596  _coords[1] = vx(1);
597  _coords[2] = vy(0);
598  _coords[3] = vy(1);
599 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
template<typename T2 >
libMesh::TypeTensor< T >::TypeTensor ( const TypeVector< T2 > &  vx,
const TypeVector< T2 > &  vy,
const TypeVector< T2 > &  vz 
)
protected

Definition at line 603 of file type_tensor.h.

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

606 {
607  libmesh_assert_equal_to (LIBMESH_DIM, 3);
608  _coords[0] = vx(0);
609  _coords[1] = vx(1);
610  _coords[2] = vx(2);
611  _coords[3] = vy(0);
612  _coords[4] = vy(1);
613  _coords[5] = vy(2);
614  _coords[6] = vz(0);
615  _coords[7] = vz(1);
616  _coords[8] = vz(2);
617 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
template<typename T2 >
libMesh::TypeTensor< T >::TypeTensor ( const TypeTensor< T2 > &  p)
inline

Copy-constructor.

Definition at line 573 of file type_tensor.h.

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

574 {
575  // copy the nodes from vector p to me
576  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
577  _coords[i] = p._coords[i];
578 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
libMesh::TypeTensor< T >::~TypeTensor ( )
inline

Destructor.

Definition at line 624 of file type_tensor.h.

625 {
626 }

Member Function Documentation

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

Add to this tensor without creating a temporary.

Definition at line 762 of file type_tensor.h.

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

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

763 {
764  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
765  _coords[i] += p._coords[i];
766 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
template<typename T2 >
void libMesh::TypeTensor< T >::add_scaled ( const TypeTensor< T2 > &  p,
const T  factor 
)
inline

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

Definition at line 773 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::TypeTensor< T >::operator=(), libMesh::System::point_hessian(), and libMesh::HPCoarsenTest::select_refinement().

774 {
775  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
776  _coords[i] += factor*p._coords[i];
777 
778 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
template<typename T2 >
void libMesh::TypeTensor< T >::assign ( const TypeTensor< T2 > &  p)
inline

Assign to this tensor without creating a temporary.

Definition at line 633 of file type_tensor.h.

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

634 {
635  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
636  _coords[i] = p._coords[i];
637 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
template<typename T2 >
CompareTypes< T, T2 >::supertype libMesh::TypeTensor< T >::contract ( const TypeTensor< T2 > &  t) const
inline

Multiply 2 tensors together, i.e. dyadic product, $ \sum_{ij} A_{ij} B_{ij} $ The tensors may contain different numeric types.

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

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

Definition at line 1171 of file type_tensor.h.

References libMesh::TypeTensor< T >::_coords, and libMesh::Parallel::sum().

Referenced by libMesh::HPCoarsenTest::add_projection(), libMesh::TensorTools::inner_product(), libMesh::TypeTensor< T >::operator=(), and libMesh::HPCoarsenTest::select_refinement().

1172 {
1173  typename CompareTypes<T,T2>::supertype sum = 0.;
1174  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1175  sum += _coords[i]*t._coords[i];
1176  return sum;
1177 }
void sum(T &r, const Communicator &comm=Communicator_World)
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
T libMesh::TypeTensor< T >::det ( ) const
inline
Returns
The determinant of the tensor.

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

Definition at line 1202 of file type_tensor.h.

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

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

1203 {
1204 #if LIBMESH_DIM == 1
1205  return _coords[0];
1206 #endif
1207 
1208 #if LIBMESH_DIM == 2
1209  return (_coords[0] * _coords[3]
1210  - _coords[1] * _coords[2]);
1211 #endif
1212 
1213 #if LIBMESH_DIM == 3
1214  return (_coords[0] * _coords[4] * _coords[8]
1215  + _coords[1] * _coords[5] * _coords[6]
1216  + _coords[2] * _coords[3] * _coords[7]
1217  - _coords[0] * _coords[5] * _coords[7]
1218  - _coords[1] * _coords[3] * _coords[8]
1219  - _coords[2] * _coords[4] * _coords[6]);
1220 #endif
1221 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
TypeTensor< T > libMesh::TypeTensor< T >::inverse ( ) const
inline
Returns
The inverse of this tensor as an independent object.

Definition at line 1018 of file type_tensor.h.

References libMesh::TypeTensor< T >::_coords, A, and libMesh::TypeTensor< T >::TypeTensor().

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

1019 {
1020 #if LIBMESH_DIM == 1
1021  if (_coords[0] == static_cast<T>(0.))
1022  libmesh_convergence_failure();
1023  return TypeTensor(1. / _coords[0]);
1024 #endif
1025 
1026 #if LIBMESH_DIM == 2
1027  // Get convenient reference to this.
1028  const TypeTensor<T> & A = *this;
1029 
1030  // Use temporary variables, avoid multiple accesses
1031  T a = A(0,0), b = A(0,1),
1032  c = A(1,0), d = A(1,1);
1033 
1034  // Make sure det = ad - bc is not zero
1035  T my_det = a*d - b*c;
1036 
1037  if (my_det == static_cast<T>(0.))
1038  libmesh_convergence_failure();
1039 
1040  return TypeTensor(d/my_det, -b/my_det, -c/my_det, a/my_det);
1041 #endif
1042 
1043 #if LIBMESH_DIM == 3
1044  // Get convenient reference to this.
1045  const TypeTensor<T> & A = *this;
1046 
1047  T a11 = A(0,0), a12 = A(0,1), a13 = A(0,2),
1048  a21 = A(1,0), a22 = A(1,1), a23 = A(1,2),
1049  a31 = A(2,0), a32 = A(2,1), a33 = A(2,2);
1050 
1051  T my_det = a11*(a33*a22-a32*a23) - a21*(a33*a12-a32*a13) + a31*(a23*a12-a22*a13);
1052 
1053  if (my_det == static_cast<T>(0.))
1054  libmesh_convergence_failure();
1055 
1056  // Inline comment characters are for lining up columns.
1057  return TypeTensor( (a33*a22-a32*a23)/my_det, -(a33*a12-a32*a13)/my_det, (a23*a12-a22*a13)/my_det,
1058  -(a33*a21-a31*a23)/my_det, (a33*a11-a31*a13)/my_det, -(a23*a11-a21*a13)/my_det,
1059  (a32*a21-a31*a22)/my_det, -(a32*a11-a31*a12)/my_det, (a22*a11-a21*a12)/my_det);
1060 #endif
1061 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
static PetscErrorCode Mat * A
template<typename T >
Real libMesh::TypeTensor< T >::norm ( ) const
inline
Returns
The Frobenius norm of the tensor, i.e. the square-root of the sum of the elements squared.

Definition at line 1193 of file type_tensor.h.

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

Referenced by libMesh::System::calculate_norm(), libMesh::TypeTensor< T >::operator=(), and libMesh::TypeTensor< T >::size().

1194 {
1195  return std::sqrt(this->norm_sq());
1196 }
Real norm_sq() const
Definition: type_tensor.h:1262
template<typename T >
Real libMesh::TypeTensor< T >::norm_sq ( ) const
inline
Returns
The Frobenius norm of the tensor squared, i.e. sum of the element magnitudes squared.

Definition at line 1262 of file type_tensor.h.

References libMesh::TypeTensor< T >::_coords, libMesh::TensorTools::norm_sq(), libMesh::Real, and libMesh::Parallel::sum().

Referenced by libMesh::UniformRefinementEstimator::_estimate_error(), libMesh::System::calculate_norm(), libMesh::ExactErrorEstimator::find_squared_element_error(), libMesh::TypeTensor< T >::norm(), libMesh::TypeTensor< T >::operator=(), libMesh::HPCoarsenTest::select_refinement(), and libMesh::TypeTensor< T >::size_sq().

1263 {
1264  Real sum = 0.;
1265  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1266  sum += TensorTools::norm_sq(_coords[i]);
1267  return sum;
1268 }
void sum(T &r, const Communicator &comm=Communicator_World)
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
template<typename T >
const T & libMesh::TypeTensor< T >::operator() ( const unsigned int  i,
const unsigned int  j 
) const
inline
Returns
A const reference to the (i,j) entry of the tensor.

Definition at line 643 of file type_tensor.h.

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

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

645 {
646  libmesh_assert_less (i, 3);
647  libmesh_assert_less (j, 3);
648 
649 #if LIBMESH_DIM < 3
650  const static T my_zero = 0;
651  if (i >= LIBMESH_DIM || j >= LIBMESH_DIM)
652  return my_zero;
653 #endif
654 
655  return _coords[i*LIBMESH_DIM+j];
656 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
T & libMesh::TypeTensor< T >::operator() ( const unsigned int  i,
const unsigned int  j 
)
inline
Returns
A writable reference to the (i,j) entry of the tensor.

Definition at line 662 of file type_tensor.h.

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

664 {
665 #if LIBMESH_DIM < 3
666 
667  if (i >= LIBMESH_DIM || j >= LIBMESH_DIM)
668  libmesh_error_msg("ERROR: You are assigning to a tensor component that is out of range for the compiled LIBMESH_DIM!");
669 
670 #endif
671 
672  libmesh_assert_less (i, LIBMESH_DIM);
673  libmesh_assert_less (j, LIBMESH_DIM);
674 
675  return _coords[i*LIBMESH_DIM+j];
676 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
template<typename Scalar >
boostcopy::enable_if_c< ScalarTraits< Scalar >::value, TypeTensor< typename CompareTypes< T, Scalar >::supertype > >::type libMesh::TypeTensor< T >::operator* ( const Scalar  factor) const
inline

Multiply this tensor by a scalar value.

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

Definition at line 889 of file type_tensor.h.

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

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

890 {
891  typedef typename CompareTypes<T, Scalar>::supertype TS;
892 
893 
894 #if LIBMESH_DIM == 1
895  return TypeTensor<TS>(_coords[0]*factor);
896 #endif
897 
898 #if LIBMESH_DIM == 2
899  return TypeTensor<TS>(_coords[0]*factor,
900  _coords[1]*factor,
901  _coords[2]*factor,
902  _coords[3]*factor);
903 #endif
904 
905 #if LIBMESH_DIM == 3
906  return TypeTensor<TS>(_coords[0]*factor,
907  _coords[1]*factor,
908  _coords[2]*factor,
909  _coords[3]*factor,
910  _coords[4]*factor,
911  _coords[5]*factor,
912  _coords[6]*factor,
913  _coords[7]*factor,
914  _coords[8]*factor);
915 #endif
916 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
template<typename T2 >
TypeTensor< T > libMesh::TypeTensor< T >::operator* ( const TypeTensor< T2 > &  p) const
inline

Multiply 2 tensors together, i.e. matrix-matrix product. The tensors may contain different numeric types.

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

Definition at line 1150 of file type_tensor.h.

1151 {
1152  TypeTensor<T> returnval;
1153  for (unsigned int i=0; i<LIBMESH_DIM; i++)
1154  for (unsigned int j=0; j<LIBMESH_DIM; j++)
1155  for (unsigned int k=0; k<LIBMESH_DIM; k++)
1156  returnval(i,j) += (*this)(i,k)*p(k,j);
1157 
1158  return returnval;
1159 }
template<typename T >
template<typename T2 >
TypeVector< typename CompareTypes< T, T2 >::supertype > libMesh::TypeTensor< T >::operator* ( const TypeVector< T2 > &  p) const
inline

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

1136 {
1137  TypeVector<typename CompareTypes<T,T2>::supertype> returnval;
1138  for (unsigned int i=0; i<LIBMESH_DIM; i++)
1139  for (unsigned int j=0; j<LIBMESH_DIM; j++)
1140  returnval(i) += (*this)(i,j)*p(j);
1141 
1142  return returnval;
1143 }
template<typename T >
template<typename Scalar >
const TypeTensor< T > & libMesh::TypeTensor< T >::operator*= ( const Scalar  factor)
inline

Multiply this tensor by a scalar value in place.

Returns
A reference to *this.

Definition at line 936 of file type_tensor.h.

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

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

937 {
938  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
939  _coords[i] *= factor;
940 
941  return *this;
942 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
template<typename T2 >
TypeTensor< typename CompareTypes< T, T2 >::supertype > libMesh::TypeTensor< T >::operator+ ( const TypeTensor< T2 > &  p) const
inline

Add another tensor to this tensor.

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

Definition at line 716 of file type_tensor.h.

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

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

717 {
718 
719 #if LIBMESH_DIM == 1
720  return TypeTensor(_coords[0] + p._coords[0]);
721 #endif
722 
723 #if LIBMESH_DIM == 2
724  return TypeTensor(_coords[0] + p._coords[0],
725  _coords[1] + p._coords[1],
726  0.,
727  _coords[2] + p._coords[2],
728  _coords[3] + p._coords[3]);
729 #endif
730 
731 #if LIBMESH_DIM == 3
732  return TypeTensor(_coords[0] + p._coords[0],
733  _coords[1] + p._coords[1],
734  _coords[2] + p._coords[2],
735  _coords[3] + p._coords[3],
736  _coords[4] + p._coords[4],
737  _coords[5] + p._coords[5],
738  _coords[6] + p._coords[6],
739  _coords[7] + p._coords[7],
740  _coords[8] + p._coords[8]);
741 #endif
742 
743 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
template<typename T2 >
const TypeTensor< T > & libMesh::TypeTensor< T >::operator+= ( const TypeTensor< T2 > &  p)
inline

Add to this tensor.

Returns
A reference to *this.

Definition at line 750 of file type_tensor.h.

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

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

751 {
752  this->add (p);
753 
754  return *this;
755 }
void add(const TypeTensor< T2 > &)
Definition: type_tensor.h:762
template<typename T >
template<typename T2 >
TypeTensor< typename CompareTypes< T, T2 >::supertype > libMesh::TypeTensor< T >::operator- ( const TypeTensor< T2 > &  p) const
inline

Subtract a tensor from this tensor.

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

Definition at line 786 of file type_tensor.h.

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

787 {
788 
789 #if LIBMESH_DIM == 1
790  return TypeTensor(_coords[0] - p._coords[0]);
791 #endif
792 
793 #if LIBMESH_DIM == 2
794  return TypeTensor(_coords[0] - p._coords[0],
795  _coords[1] - p._coords[1],
796  0.,
797  _coords[2] - p._coords[2],
798  _coords[3] - p._coords[3]);
799 #endif
800 
801 #if LIBMESH_DIM == 3
802  return TypeTensor(_coords[0] - p._coords[0],
803  _coords[1] - p._coords[1],
804  _coords[2] - p._coords[2],
805  _coords[3] - p._coords[3],
806  _coords[4] - p._coords[4],
807  _coords[5] - p._coords[5],
808  _coords[6] - p._coords[6],
809  _coords[7] - p._coords[7],
810  _coords[8] - p._coords[8]);
811 #endif
812 
813 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
TypeTensor< T > libMesh::TypeTensor< T >::operator- ( ) const
inline
Returns
The negative of this tensor in a separate copy.

Definition at line 853 of file type_tensor.h.

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

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

854 {
855 
856 #if LIBMESH_DIM == 1
857  return TypeTensor(-_coords[0]);
858 #endif
859 
860 #if LIBMESH_DIM == 2
861  return TypeTensor(-_coords[0],
862  -_coords[1],
863  -_coords[2],
864  -_coords[3]);
865 #endif
866 
867 #if LIBMESH_DIM == 3
868  return TypeTensor(-_coords[0],
869  -_coords[1],
870  -_coords[2],
871  -_coords[3],
872  -_coords[4],
873  -_coords[5],
874  -_coords[6],
875  -_coords[7],
876  -_coords[8]);
877 #endif
878 
879 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
template<typename T2 >
const TypeTensor< T > & libMesh::TypeTensor< T >::operator-= ( const TypeTensor< T2 > &  p)
inline

Subtract from this tensor.

Returns
A reference to *this.

Definition at line 820 of file type_tensor.h.

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

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

821 {
822  this->subtract (p);
823 
824  return *this;
825 }
void subtract(const TypeTensor< T2 > &)
Definition: type_tensor.h:832
template<typename T >
template<typename Scalar >
boostcopy::enable_if_c< ScalarTraits< Scalar >::value, TypeTensor< typename CompareTypes< T, Scalar >::supertype > >::type libMesh::TypeTensor< T >::operator/ ( const Scalar  factor) const
inline

Divide each entry of this tensor by a scalar value.

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

Definition at line 953 of file type_tensor.h.

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

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

954 {
955  libmesh_assert_not_equal_to (factor, static_cast<T>(0.));
956 
957  typedef typename CompareTypes<T, Scalar>::supertype TS;
958 
959 #if LIBMESH_DIM == 1
960  return TypeTensor<TS>(_coords[0]/factor);
961 #endif
962 
963 #if LIBMESH_DIM == 2
964  return TypeTensor<TS>(_coords[0]/factor,
965  _coords[1]/factor,
966  _coords[2]/factor,
967  _coords[3]/factor);
968 #endif
969 
970 #if LIBMESH_DIM == 3
971  return TypeTensor<TS>(_coords[0]/factor,
972  _coords[1]/factor,
973  _coords[2]/factor,
974  _coords[3]/factor,
975  _coords[4]/factor,
976  _coords[5]/factor,
977  _coords[6]/factor,
978  _coords[7]/factor,
979  _coords[8]/factor);
980 #endif
981 
982 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
const TypeTensor< T > & libMesh::TypeTensor< T >::operator/= ( const T  factor)
inline

Divide each entry of this tensor by a scalar value.

Returns
A reference to *this.

Definition at line 1118 of file type_tensor.h.

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

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

1119 {
1120  libmesh_assert_not_equal_to (factor, static_cast<T>(0.));
1121 
1122  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1123  _coords[i] /= factor;
1124 
1125  return *this;
1126 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<>
bool libMesh::TypeTensor< Real >::operator< ( const TypeTensor< Real > &  rhs) const

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 }
template<>
bool libMesh::TypeTensor< Complex >::operator< ( const TypeTensor< Complex > &  rhs) const

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 }
template<typename T>
bool libMesh::TypeTensor< T >::operator< ( const TypeTensor< T > &  rhs) const
Returns
true if this tensor is "less" than another. Useful for sorting.
template<typename T>
template<typename Scalar >
boostcopy::enable_if_c< ScalarTraits<Scalar>::value, TypeTensor &>::type libMesh::TypeTensor< T >::operator= ( const Scalar &  libmesh_dbg_varp)
inline
template<typename T >
bool libMesh::TypeTensor< T >::operator== ( const TypeTensor< T > &  rhs) const
inline
Returns
true if two tensors are equal, false otherwise.

Definition at line 1274 of file type_tensor.h.

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

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

1275 {
1276 #if LIBMESH_DIM == 1
1277  return (std::abs(_coords[0] - rhs._coords[0])
1278  < TOLERANCE);
1279 #endif
1280 
1281 #if LIBMESH_DIM == 2
1282  return ((std::abs(_coords[0] - rhs._coords[0]) +
1283  std::abs(_coords[1] - rhs._coords[1]) +
1284  std::abs(_coords[2] - rhs._coords[2]) +
1285  std::abs(_coords[3] - rhs._coords[3]))
1286  < 4.*TOLERANCE);
1287 #endif
1288 
1289 #if LIBMESH_DIM == 3
1290  return ((std::abs(_coords[0] - rhs._coords[0]) +
1291  std::abs(_coords[1] - rhs._coords[1]) +
1292  std::abs(_coords[2] - rhs._coords[2]) +
1293  std::abs(_coords[3] - rhs._coords[3]) +
1294  std::abs(_coords[4] - rhs._coords[4]) +
1295  std::abs(_coords[5] - rhs._coords[5]) +
1296  std::abs(_coords[6] - rhs._coords[6]) +
1297  std::abs(_coords[7] - rhs._coords[7]) +
1298  std::abs(_coords[8] - rhs._coords[8]))
1299  < 9.*TOLERANCE);
1300 #endif
1301 
1302 }
double abs(double a)
static const Real TOLERANCE
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<>
bool libMesh::TypeTensor< Real >::operator> ( const TypeTensor< Real > &  rhs) const

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 }
template<>
bool libMesh::TypeTensor< Complex >::operator> ( const TypeTensor< Complex > &  rhs) const

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 }
template<typename T>
bool libMesh::TypeTensor< T >::operator> ( const TypeTensor< T > &  rhs) const
Returns
true if this tensor is "greater" than another.

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

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

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

Definition at line 39 of file type_tensor.C.

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

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 }
template<typename T >
TypeVector< T > libMesh::TypeTensor< T >::row ( const unsigned int  r) const
inline
Returns
A copy of one row of the tensor as a TypeVector.

Definition at line 702 of file type_tensor.h.

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

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

703 {
704  TypeVector<T> return_vector;
705 
706  for (unsigned int j=0; j<LIBMESH_DIM; j++)
707  return_vector._coords[j] = _coords[r*LIBMESH_DIM + j];
708 
709  return return_vector;
710 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
Real libMesh::TypeTensor< T >::size ( ) const
inline
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 1183 of file type_tensor.h.

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

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

1184 {
1185  libmesh_deprecated();
1186  return this->norm();
1187 }
Real norm() const
Definition: type_tensor.h:1193
template<typename T >
Real libMesh::TypeTensor< T >::size_sq ( ) const
inline
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 1252 of file type_tensor.h.

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

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

1253 {
1254  libmesh_deprecated();
1255  return this->norm_sq();
1256 }
Real norm_sq() const
Definition: type_tensor.h:1262
template<typename T >
ConstTypeTensorColumn< T > libMesh::TypeTensor< T >::slice ( const unsigned int  i) const
inline
Returns
A proxy for the $ i^{th} $ column of the tensor.

Definition at line 682 of file type_tensor.h.

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

683 {
684  libmesh_assert_less (i, LIBMESH_DIM);
685  return ConstTypeTensorColumn<T>(*this, i);
686 }
template<typename T >
TypeTensorColumn< T > libMesh::TypeTensor< T >::slice ( const unsigned int  i)
inline
Returns
A writable proxy for the $ i^{th} $ column of the tensor.

Definition at line 692 of file type_tensor.h.

693 {
694  libmesh_assert_less (i, LIBMESH_DIM);
695  return TypeTensorColumn<T>(*this, i);
696 }
template<typename T >
void libMesh::TypeTensor< T >::solve ( const TypeVector< T > &  b,
TypeVector< T > &  x 
) const
inline

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

Returns
The solution in the x vector.

Definition at line 1067 of file type_tensor.h.

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

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

1068 {
1069 #if LIBMESH_DIM == 1
1070  if (_coords[0] == static_cast<T>(0.))
1071  libmesh_convergence_failure();
1072  x(0) = b(0) / _coords[0];
1073 #endif
1074 
1075 #if LIBMESH_DIM == 2
1076  T my_det = _coords[0]*_coords[3] - _coords[1]*_coords[2];
1077 
1078  if (my_det == static_cast<T>(0.))
1079  libmesh_convergence_failure();
1080 
1081  T my_det_inv = 1./my_det;
1082 
1083  x(0) = my_det_inv*( _coords[3]*b(0) - _coords[1]*b(1));
1084  x(1) = my_det_inv*(-_coords[2]*b(0) + _coords[0]*b(1));
1085 #endif
1086 
1087 #if LIBMESH_DIM == 3
1088  T my_det =
1089  // a11*(a33 *a22 - a32 *a23)
1090  _coords[0]*(_coords[8]*_coords[4] - _coords[7]*_coords[5])
1091  // -a21*(a33 *a12 - a32 *a13)
1092  -_coords[3]*(_coords[8]*_coords[1] - _coords[7]*_coords[2]) +
1093  // +a31*(a23 *a12 - a22 *a13)
1094  +_coords[6]*(_coords[5]*_coords[1] - _coords[4]*_coords[2]);
1095 
1096  if (my_det == static_cast<T>(0.))
1097  libmesh_convergence_failure();
1098 
1099  T my_det_inv = 1./my_det;
1100  x(0) = my_det_inv*((_coords[8]*_coords[4] - _coords[7]*_coords[5])*b(0) -
1101  (_coords[8]*_coords[1] - _coords[7]*_coords[2])*b(1) +
1102  (_coords[5]*_coords[1] - _coords[4]*_coords[2])*b(2));
1103 
1104  x(1) = my_det_inv*((_coords[6]*_coords[5] - _coords[8]*_coords[3])*b(0) +
1105  (_coords[8]*_coords[0] - _coords[6]*_coords[2])*b(1) -
1106  (_coords[5]*_coords[0] - _coords[3]*_coords[2])*b(2));
1107 
1108  x(2) = my_det_inv*((_coords[7]*_coords[3] - _coords[6]*_coords[4])*b(0) -
1109  (_coords[7]*_coords[0] - _coords[6]*_coords[1])*b(1) +
1110  (_coords[4]*_coords[0] - _coords[3]*_coords[1])*b(2));
1111 #endif
1112 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
PetscErrorCode Vec x
template<typename T >
template<typename T2 >
void libMesh::TypeTensor< T >::subtract ( const TypeTensor< T2 > &  p)
inline

Subtract from this tensor without creating a temporary.

Definition at line 832 of file type_tensor.h.

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

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

833 {
834  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
835  _coords[i] -= p._coords[i];
836 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
template<typename T2 >
void libMesh::TypeTensor< T >::subtract_scaled ( const TypeTensor< T2 > &  p,
const T  factor 
)
inline

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

Definition at line 843 of file type_tensor.h.

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

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

844 {
845  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
846  _coords[i] -= factor*p._coords[i];
847 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
T libMesh::TypeTensor< T >::tr ( ) const
inline
Returns
The trace of the tensor.

Definition at line 1225 of file type_tensor.h.

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

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

1226 {
1227 #if LIBMESH_DIM == 1
1228  return _coords[0];
1229 #endif
1230 
1231 #if LIBMESH_DIM == 2
1232  return _coords[0] + _coords[3];
1233 #endif
1234 
1235 #if LIBMESH_DIM == 3
1236  return _coords[0] + _coords[4] + _coords[8];
1237 #endif
1238 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
TypeTensor< T > libMesh::TypeTensor< T >::transpose ( ) const
inline
Returns
The transpose of this tensor (with complex numbers not conjugated).

Definition at line 988 of file type_tensor.h.

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

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

989 {
990 #if LIBMESH_DIM == 1
991  return TypeTensor(_coords[0]);
992 #endif
993 
994 #if LIBMESH_DIM == 2
995  return TypeTensor(_coords[0],
996  _coords[2],
997  _coords[1],
998  _coords[3]);
999 #endif
1000 
1001 #if LIBMESH_DIM == 3
1002  return TypeTensor(_coords[0],
1003  _coords[3],
1004  _coords[6],
1005  _coords[1],
1006  _coords[4],
1007  _coords[7],
1008  _coords[2],
1009  _coords[5],
1010  _coords[8]);
1011 #endif
1012 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413
template<typename T >
void libMesh::TypeTensor< T >::write_unformatted ( std::ostream &  out,
const bool  newline = true 
) const

Unformatted print to the stream out. Simply prints the elements of the tensor separated by spaces and newlines.

Definition at line 83 of file type_tensor.C.

References libMesh::libmesh_assert().

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 }
libmesh_assert(j)
template<typename T >
void libMesh::TypeTensor< T >::zero ( )
inline

Set all entries of the tensor to 0.

Definition at line 1242 of file type_tensor.h.

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

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

1243 {
1244  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1245  _coords[i] = 0.;
1246 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:413

Friends And Related Function Documentation

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

Formatted print as above but supports the syntax:

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

Definition at line 396 of file type_tensor.h.

397  {
398  t.print(os);
399  return os;
400  }
template<typename T>
template<typename T2 >
friend class TypeTensor
friend

Definition at line 52 of file type_tensor.h.

Member Data Documentation


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