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
 

Public Attributes

_coords [LIBMESH_DIM *LIBMESH_DIM]
 

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)
 

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 450 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().

451 {
452  _coords[0] = 0;
453 
454 #if LIBMESH_DIM > 1
455  _coords[1] = 0;
456  _coords[2] = 0;
457  _coords[3] = 0;
458 #endif
459 
460 #if LIBMESH_DIM > 2
461  _coords[4] = 0;
462  _coords[5] = 0;
463  _coords[6] = 0;
464  _coords[7] = 0;
465  _coords[8] = 0;
466 #endif
467 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
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 473 of file type_tensor.h.

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

482 {
483  _coords[0] = xx;
484 
485 #if LIBMESH_DIM == 2
486  _coords[1] = xy;
487  _coords[2] = yx;
488  _coords[3] = yy;
489 #endif
490 
491 #if LIBMESH_DIM == 3
492  _coords[1] = xy;
493  _coords[2] = xz;
494  _coords[3] = yx;
495  _coords[4] = yy;
496  _coords[5] = yz;
497  _coords[6] = zx;
498  _coords[7] = zy;
499  _coords[8] = zz;
500 #endif
501 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
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 507 of file type_tensor.h.

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

518 {
519  _coords[0] = xx;
520 
521 #if LIBMESH_DIM == 2
522  _coords[1] = xy;
523  _coords[2] = yx;
524  _coords[3] = yy;
525 #endif
526 
527 #if LIBMESH_DIM == 3
528  _coords[1] = xy;
529  _coords[2] = xz;
530  _coords[3] = yx;
531  _coords[4] = yy;
532  _coords[5] = yz;
533  _coords[6] = zx;
534  _coords[7] = zy;
535  _coords[8] = zz;
536 #endif
537 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
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 554 of file type_tensor.h.

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

555 {
556  libmesh_assert_equal_to (LIBMESH_DIM, 1);
557  _coords[0] = vx(0);
558 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
template<typename T >
template<typename T2 >
libMesh::TypeTensor< T >::TypeTensor ( const TypeVector< T2 > &  vx,
const TypeVector< T2 > &  vy 
)
protected

Definition at line 562 of file type_tensor.h.

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

564 {
565  libmesh_assert_equal_to (LIBMESH_DIM, 2);
566  _coords[0] = vx(0);
567  _coords[1] = vx(1);
568  _coords[2] = vy(0);
569  _coords[3] = vy(1);
570 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
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 574 of file type_tensor.h.

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

577 {
578  libmesh_assert_equal_to (LIBMESH_DIM, 3);
579  _coords[0] = vx(0);
580  _coords[1] = vx(1);
581  _coords[2] = vx(2);
582  _coords[3] = vy(0);
583  _coords[4] = vy(1);
584  _coords[5] = vy(2);
585  _coords[6] = vz(0);
586  _coords[7] = vz(1);
587  _coords[8] = vz(2);
588 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
template<typename T >
template<typename T2 >
libMesh::TypeTensor< T >::TypeTensor ( const TypeTensor< T2 > &  p)
inline

Copy-constructor.

Definition at line 544 of file type_tensor.h.

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

545 {
546  // copy the nodes from vector p to me
547  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
548  _coords[i] = p._coords[i];
549 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
template<typename T >
libMesh::TypeTensor< T >::~TypeTensor ( )
inline

Destructor.

Definition at line 595 of file type_tensor.h.

596 {
597 }

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

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

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

734 {
735  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
736  _coords[i] += p._coords[i];
737 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
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 744 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().

745 {
746  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
747  _coords[i] += factor*p._coords[i];
748 
749 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
template<typename T >
template<typename T2 >
void libMesh::TypeTensor< T >::assign ( const TypeTensor< T2 > &  p)
inline

Assign to a tensor without creating a temporary.

Definition at line 604 of file type_tensor.h.

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

605 {
606  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
607  _coords[i] = p._coords[i];
608 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
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 Aij*Bij. The tensors may be of different types.

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

Definition at line 1142 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().

1143 {
1144  typename CompareTypes<T,T2>::supertype sum = 0.;
1145  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1146  sum += _coords[i]*t._coords[i];
1147  return sum;
1148 }
void sum(T &r, const Communicator &comm=Communicator_World)
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
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 1173 of file type_tensor.h.

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

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

1174 {
1175 #if LIBMESH_DIM == 1
1176  return _coords[0];
1177 #endif
1178 
1179 #if LIBMESH_DIM == 2
1180  return (_coords[0] * _coords[3]
1181  - _coords[1] * _coords[2]);
1182 #endif
1183 
1184 #if LIBMESH_DIM == 3
1185  return (_coords[0] * _coords[4] * _coords[8]
1186  + _coords[1] * _coords[5] * _coords[6]
1187  + _coords[2] * _coords[3] * _coords[7]
1188  - _coords[0] * _coords[5] * _coords[7]
1189  - _coords[1] * _coords[3] * _coords[8]
1190  - _coords[2] * _coords[4] * _coords[6]);
1191 #endif
1192 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
template<typename T >
TypeTensor< T > libMesh::TypeTensor< T >::inverse ( ) const
inline

Returns the inverse of the tensor.

Definition at line 989 of file type_tensor.h.

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

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

990 {
991 #if LIBMESH_DIM == 1
992  if (_coords[0] == static_cast<T>(0.))
993  libmesh_convergence_failure();
994  return TypeTensor(1. / _coords[0]);
995 #endif
996 
997 #if LIBMESH_DIM == 2
998  // Get convenient reference to this.
999  const TypeTensor<T> & A = *this;
1000 
1001  // Use temporary variables, avoid multiple accesses
1002  T a = A(0,0), b = A(0,1),
1003  c = A(1,0), d = A(1,1);
1004 
1005  // Make sure det = ad - bc is not zero
1006  T my_det = a*d - b*c;
1007 
1008  if (my_det == static_cast<T>(0.))
1009  libmesh_convergence_failure();
1010 
1011  return TypeTensor(d/my_det, -b/my_det, -c/my_det, a/my_det);
1012 #endif
1013 
1014 #if LIBMESH_DIM == 3
1015  // Get convenient reference to this.
1016  const TypeTensor<T> & A = *this;
1017 
1018  T a11 = A(0,0), a12 = A(0,1), a13 = A(0,2),
1019  a21 = A(1,0), a22 = A(1,1), a23 = A(1,2),
1020  a31 = A(2,0), a32 = A(2,1), a33 = A(2,2);
1021 
1022  T my_det = a11*(a33*a22-a32*a23) - a21*(a33*a12-a32*a13) + a31*(a23*a12-a22*a13);
1023 
1024  if (my_det == static_cast<T>(0.))
1025  libmesh_convergence_failure();
1026 
1027  // Inline comment characters are for lining up columns.
1028  return TypeTensor( (a33*a22-a32*a23)/my_det, -(a33*a12-a32*a13)/my_det, (a23*a12-a22*a13)/my_det,
1029  -(a33*a21-a31*a23)/my_det, (a33*a11-a31*a13)/my_det, -(a23*a11-a21*a13)/my_det,
1030  (a32*a21-a31*a22)/my_det, -(a32*a11-a31*a12)/my_det, (a22*a11-a21*a12)/my_det);
1031 #endif
1032 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
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 1164 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().

1165 {
1166  return std::sqrt(this->norm_sq());
1167 }
Real norm_sq() const
Definition: type_tensor.h:1233
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 1233 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().

1234 {
1235  Real sum = 0.;
1236  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1237  sum += TensorTools::norm_sq(_coords[i]);
1238  return sum;
1239 }
void sum(T &r, const Communicator &comm=Communicator_World)
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
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

Return the $ i,j^{th} $ element of the tensor.

Definition at line 614 of file type_tensor.h.

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

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

616 {
617  libmesh_assert_less (i, 3);
618  libmesh_assert_less (j, 3);
619 
620 #if LIBMESH_DIM < 3
621  const static T my_zero = 0;
622  if (i >= LIBMESH_DIM || j >= LIBMESH_DIM)
623  return my_zero;
624 #endif
625 
626  return _coords[i*LIBMESH_DIM+j];
627 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
template<typename T >
T & libMesh::TypeTensor< T >::operator() ( const unsigned int  i,
const unsigned int  j 
)
inline

Return a writeable reference to the $ i,j^{th} $ element of the tensor.

Definition at line 633 of file type_tensor.h.

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

635 {
636 #if LIBMESH_DIM < 3
637 
638  if (i >= LIBMESH_DIM || j >= LIBMESH_DIM)
639  libmesh_error_msg("ERROR: You are assigning to a tensor component that is out of range for the compiled LIBMESH_DIM!");
640 
641 #endif
642 
643  libmesh_assert_less (i, LIBMESH_DIM);
644  libmesh_assert_less (j, LIBMESH_DIM);
645 
646  return _coords[i*LIBMESH_DIM+j];
647 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
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 a tensor by a number, i.e. scale.

Definition at line 860 of file type_tensor.h.

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

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

861 {
862  typedef typename CompareTypes<T, Scalar>::supertype TS;
863 
864 
865 #if LIBMESH_DIM == 1
866  return TypeTensor<TS>(_coords[0]*factor);
867 #endif
868 
869 #if LIBMESH_DIM == 2
870  return TypeTensor<TS>(_coords[0]*factor,
871  _coords[1]*factor,
872  _coords[2]*factor,
873  _coords[3]*factor);
874 #endif
875 
876 #if LIBMESH_DIM == 3
877  return TypeTensor<TS>(_coords[0]*factor,
878  _coords[1]*factor,
879  _coords[2]*factor,
880  _coords[3]*factor,
881  _coords[4]*factor,
882  _coords[5]*factor,
883  _coords[6]*factor,
884  _coords[7]*factor,
885  _coords[8]*factor);
886 #endif
887 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
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 product. The tensors may be of different types.

Definition at line 1121 of file type_tensor.h.

1122 {
1123  TypeTensor<T> returnval;
1124  for (unsigned int i=0; i<LIBMESH_DIM; i++)
1125  for (unsigned int j=0; j<LIBMESH_DIM; j++)
1126  for (unsigned int k=0; k<LIBMESH_DIM; k++)
1127  returnval(i,j) += (*this)(i,k)*p(k,j);
1128 
1129  return returnval;
1130 }
template<typename T >
template<typename T2 >
TypeVector< typename CompareTypes< T, T2 >::supertype > libMesh::TypeTensor< T >::operator* ( const TypeVector< T2 > &  p) const
inline

Multiply a tensor and vector together, i.e. matrix-vector product. The tensor and vector may be of different types.

Definition at line 1106 of file type_tensor.h.

1107 {
1108  TypeVector<typename CompareTypes<T,T2>::supertype> returnval;
1109  for (unsigned int i=0; i<LIBMESH_DIM; i++)
1110  for (unsigned int j=0; j<LIBMESH_DIM; j++)
1111  returnval(i) += (*this)(i,j)*p(j);
1112 
1113  return returnval;
1114 }
template<typename T >
template<typename Scalar >
const TypeTensor< T > & libMesh::TypeTensor< T >::operator*= ( const Scalar  factor)
inline

Multiply this tensor by a number, i.e. scale.

Definition at line 907 of file type_tensor.h.

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

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

908 {
909  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
910  _coords[i] *= factor;
911 
912  return *this;
913 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
template<typename T >
template<typename T2 >
TypeTensor< typename CompareTypes< T, T2 >::supertype > libMesh::TypeTensor< T >::operator+ ( const TypeTensor< T2 > &  p) const
inline

Add two tensors.

Definition at line 687 of file type_tensor.h.

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

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

688 {
689 
690 #if LIBMESH_DIM == 1
691  return TypeTensor(_coords[0] + p._coords[0]);
692 #endif
693 
694 #if LIBMESH_DIM == 2
695  return TypeTensor(_coords[0] + p._coords[0],
696  _coords[1] + p._coords[1],
697  0.,
698  _coords[2] + p._coords[2],
699  _coords[3] + p._coords[3]);
700 #endif
701 
702 #if LIBMESH_DIM == 3
703  return TypeTensor(_coords[0] + p._coords[0],
704  _coords[1] + p._coords[1],
705  _coords[2] + p._coords[2],
706  _coords[3] + p._coords[3],
707  _coords[4] + p._coords[4],
708  _coords[5] + p._coords[5],
709  _coords[6] + p._coords[6],
710  _coords[7] + p._coords[7],
711  _coords[8] + p._coords[8]);
712 #endif
713 
714 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
template<typename T >
template<typename T2 >
const TypeTensor< T > & libMesh::TypeTensor< T >::operator+= ( const TypeTensor< T2 > &  p)
inline

Add to this tensor.

Definition at line 721 of file type_tensor.h.

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

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

722 {
723  this->add (p);
724 
725  return *this;
726 }
void add(const TypeTensor< T2 > &)
Definition: type_tensor.h:733
template<typename T >
template<typename T2 >
TypeTensor< typename CompareTypes< T, T2 >::supertype > libMesh::TypeTensor< T >::operator- ( const TypeTensor< T2 > &  p) const
inline

Subtract two tensors.

Definition at line 757 of file type_tensor.h.

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

758 {
759 
760 #if LIBMESH_DIM == 1
761  return TypeTensor(_coords[0] - p._coords[0]);
762 #endif
763 
764 #if LIBMESH_DIM == 2
765  return TypeTensor(_coords[0] - p._coords[0],
766  _coords[1] - p._coords[1],
767  0.,
768  _coords[2] - p._coords[2],
769  _coords[3] - p._coords[3]);
770 #endif
771 
772 #if LIBMESH_DIM == 3
773  return TypeTensor(_coords[0] - p._coords[0],
774  _coords[1] - p._coords[1],
775  _coords[2] - p._coords[2],
776  _coords[3] - p._coords[3],
777  _coords[4] - p._coords[4],
778  _coords[5] - p._coords[5],
779  _coords[6] - p._coords[6],
780  _coords[7] - p._coords[7],
781  _coords[8] - p._coords[8]);
782 #endif
783 
784 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
template<typename T >
TypeTensor< T > libMesh::TypeTensor< T >::operator- ( ) const
inline

Return the opposite of a tensor

Definition at line 824 of file type_tensor.h.

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

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

825 {
826 
827 #if LIBMESH_DIM == 1
828  return TypeTensor(-_coords[0]);
829 #endif
830 
831 #if LIBMESH_DIM == 2
832  return TypeTensor(-_coords[0],
833  -_coords[1],
834  -_coords[2],
835  -_coords[3]);
836 #endif
837 
838 #if LIBMESH_DIM == 3
839  return TypeTensor(-_coords[0],
840  -_coords[1],
841  -_coords[2],
842  -_coords[3],
843  -_coords[4],
844  -_coords[5],
845  -_coords[6],
846  -_coords[7],
847  -_coords[8]);
848 #endif
849 
850 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
template<typename T >
template<typename T2 >
const TypeTensor< T > & libMesh::TypeTensor< T >::operator-= ( const TypeTensor< T2 > &  p)
inline

Subtract from this tensor.

Definition at line 791 of file type_tensor.h.

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

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

792 {
793  this->subtract (p);
794 
795  return *this;
796 }
void subtract(const TypeTensor< T2 > &)
Definition: type_tensor.h:803
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 a tensor by a number, i.e. scale.

Definition at line 924 of file type_tensor.h.

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

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

925 {
926  libmesh_assert_not_equal_to (factor, static_cast<T>(0.));
927 
928  typedef typename CompareTypes<T, Scalar>::supertype TS;
929 
930 #if LIBMESH_DIM == 1
931  return TypeTensor<TS>(_coords[0]/factor);
932 #endif
933 
934 #if LIBMESH_DIM == 2
935  return TypeTensor<TS>(_coords[0]/factor,
936  _coords[1]/factor,
937  _coords[2]/factor,
938  _coords[3]/factor);
939 #endif
940 
941 #if LIBMESH_DIM == 3
942  return TypeTensor<TS>(_coords[0]/factor,
943  _coords[1]/factor,
944  _coords[2]/factor,
945  _coords[3]/factor,
946  _coords[4]/factor,
947  _coords[5]/factor,
948  _coords[6]/factor,
949  _coords[7]/factor,
950  _coords[8]/factor);
951 #endif
952 
953 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
template<typename T >
const TypeTensor< T > & libMesh::TypeTensor< T >::operator/= ( const T  factor)
inline

Divide this tensor by a number, i.e. scale.

Definition at line 1089 of file type_tensor.h.

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

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

1090 {
1091  libmesh_assert_not_equal_to (factor, static_cast<T>(0.));
1092 
1093  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1094  _coords[i] /= factor;
1095 
1096  return *this;
1097 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
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 valued.

Definition at line 1245 of file type_tensor.h.

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

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

1246 {
1247 #if LIBMESH_DIM == 1
1248  return (std::abs(_coords[0] - rhs._coords[0])
1249  < TOLERANCE);
1250 #endif
1251 
1252 #if LIBMESH_DIM == 2
1253  return ((std::abs(_coords[0] - rhs._coords[0]) +
1254  std::abs(_coords[1] - rhs._coords[1]) +
1255  std::abs(_coords[2] - rhs._coords[2]) +
1256  std::abs(_coords[3] - rhs._coords[3]))
1257  < 4.*TOLERANCE);
1258 #endif
1259 
1260 #if LIBMESH_DIM == 3
1261  return ((std::abs(_coords[0] - rhs._coords[0]) +
1262  std::abs(_coords[1] - rhs._coords[1]) +
1263  std::abs(_coords[2] - rhs._coords[2]) +
1264  std::abs(_coords[3] - rhs._coords[3]) +
1265  std::abs(_coords[4] - rhs._coords[4]) +
1266  std::abs(_coords[5] - rhs._coords[5]) +
1267  std::abs(_coords[6] - rhs._coords[6]) +
1268  std::abs(_coords[7] - rhs._coords[7]) +
1269  std::abs(_coords[8] - rhs._coords[8]))
1270  < 9.*TOLERANCE);
1271 #endif
1272 
1273 }
double abs(double a)
static const Real TOLERANCE
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
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, by default 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

Return one row of the tensor as a TypeVector.

Definition at line 673 of file type_tensor.h.

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

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

674 {
675  TypeVector<T> return_vector;
676 
677  for(unsigned int j=0; j<LIBMESH_DIM; j++)
678  return_vector._coords[j] = _coords[r*LIBMESH_DIM + j];
679 
680  return return_vector;
681 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
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. This function is deprecated, used norm() instead.

Definition at line 1154 of file type_tensor.h.

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

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

1155 {
1156  libmesh_deprecated();
1157  return this->norm();
1158 }
Real norm() const
Definition: type_tensor.h:1164
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. This function is deprecated, used norm() instead.

Definition at line 1223 of file type_tensor.h.

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

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

1224 {
1225  libmesh_deprecated();
1226  return this->norm_sq();
1227 }
Real norm_sq() const
Definition: type_tensor.h:1233
template<typename T >
ConstTypeTensorColumn< T > libMesh::TypeTensor< T >::slice ( const unsigned int  i) const
inline

Return a proxy for the $ i^{th} $ column of the tensor.

Definition at line 653 of file type_tensor.h.

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

654 {
655  libmesh_assert_less (i, LIBMESH_DIM);
656  return ConstTypeTensorColumn<T>(*this, i);
657 }
template<typename T >
TypeTensorColumn< T > libMesh::TypeTensor< T >::slice ( const unsigned int  i)
inline

Return a writeable proxy for the $ i^{th} $ column of the tensor.

Definition at line 663 of file type_tensor.h.

664 {
665  libmesh_assert_less (i, LIBMESH_DIM);
666  return TypeTensorColumn<T>(*this, i);
667 }
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.

Definition at line 1038 of file type_tensor.h.

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

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

1039 {
1040 #if LIBMESH_DIM == 1
1041  if (_coords[0] == static_cast<T>(0.))
1042  libmesh_convergence_failure();
1043  x(0) = b(0) / _coords[0];
1044 #endif
1045 
1046 #if LIBMESH_DIM == 2
1047  T my_det = _coords[0]*_coords[3] - _coords[1]*_coords[2];
1048 
1049  if (my_det == static_cast<T>(0.))
1050  libmesh_convergence_failure();
1051 
1052  T my_det_inv = 1./my_det;
1053 
1054  x(0) = my_det_inv*( _coords[3]*b(0) - _coords[1]*b(1));
1055  x(1) = my_det_inv*(-_coords[2]*b(0) + _coords[0]*b(1));
1056 #endif
1057 
1058 #if LIBMESH_DIM == 3
1059  T my_det =
1060  // a11*(a33 *a22 - a32 *a23)
1061  _coords[0]*(_coords[8]*_coords[4] - _coords[7]*_coords[5])
1062  // -a21*(a33 *a12 - a32 *a13)
1063  -_coords[3]*(_coords[8]*_coords[1] - _coords[7]*_coords[2]) +
1064  // +a31*(a23 *a12 - a22 *a13)
1065  +_coords[6]*(_coords[5]*_coords[1] - _coords[4]*_coords[2]);
1066 
1067  if (my_det == static_cast<T>(0.))
1068  libmesh_convergence_failure();
1069 
1070  T my_det_inv = 1./my_det;
1071  x(0) = my_det_inv*((_coords[8]*_coords[4] - _coords[7]*_coords[5])*b(0) -
1072  (_coords[8]*_coords[1] - _coords[7]*_coords[2])*b(1) +
1073  (_coords[5]*_coords[1] - _coords[4]*_coords[2])*b(2));
1074 
1075  x(1) = my_det_inv*((_coords[6]*_coords[5] - _coords[8]*_coords[3])*b(0) +
1076  (_coords[8]*_coords[0] - _coords[6]*_coords[2])*b(1) -
1077  (_coords[5]*_coords[0] - _coords[3]*_coords[2])*b(2));
1078 
1079  x(2) = my_det_inv*((_coords[7]*_coords[3] - _coords[6]*_coords[4])*b(0) -
1080  (_coords[7]*_coords[0] - _coords[6]*_coords[1])*b(1) +
1081  (_coords[4]*_coords[0] - _coords[3]*_coords[1])*b(2));
1082 #endif
1083 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
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 803 of file type_tensor.h.

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

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

804 {
805  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
806  _coords[i] -= p._coords[i];
807 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
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 814 of file type_tensor.h.

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

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

815 {
816  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
817  _coords[i] -= factor*p._coords[i];
818 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
template<typename T >
T libMesh::TypeTensor< T >::tr ( ) const
inline

Returns the trace of the tensor.

Definition at line 1196 of file type_tensor.h.

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

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

1197 {
1198 #if LIBMESH_DIM == 1
1199  return _coords[0];
1200 #endif
1201 
1202 #if LIBMESH_DIM == 2
1203  return _coords[0] + _coords[3];
1204 #endif
1205 
1206 #if LIBMESH_DIM == 3
1207  return _coords[0] + _coords[4] + _coords[8];
1208 #endif
1209 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
template<typename T >
TypeTensor< T > libMesh::TypeTensor< T >::transpose ( ) const
inline

The transpose (with complex numbers not conjugated) of the tensor.

Definition at line 959 of file type_tensor.h.

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

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

960 {
961 #if LIBMESH_DIM == 1
962  return TypeTensor(_coords[0]);
963 #endif
964 
965 #if LIBMESH_DIM == 2
966  return TypeTensor(_coords[0],
967  _coords[2],
968  _coords[1],
969  _coords[3]);
970 #endif
971 
972 #if LIBMESH_DIM == 3
973  return TypeTensor(_coords[0],
974  _coords[3],
975  _coords[6],
976  _coords[1],
977  _coords[4],
978  _coords[7],
979  _coords[2],
980  _coords[5],
981  _coords[8]);
982 #endif
983 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384
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

Zero the tensor in any dimension.

Definition at line 1213 of file type_tensor.h.

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

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

1214 {
1215  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1216  _coords[i] = 0.;
1217 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
Definition: type_tensor.h:384

Friends And Related Function Documentation

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

Formatted print as above but allows you to do std::cout << t << std::endl;

Definition at line 367 of file type_tensor.h.

368  {
369  t.print(os);
370  return os;
371  }
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: