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

This class defines a tensor in LIBMESH_DIM dimensional space of type T. More...

#include <tensor_tools.h>

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

Public Types

typedef T value_type
 Helper typedef for C++98 generic programming. More...
 

Public Member Functions

template<typename T2 >
 TypeTensor (const TypeTensor< T2 > &p)
 Copy-constructor. More...
 
 ~TypeTensor ()
 Destructor. More...
 
template<typename T2 >
void assign (const TypeTensor< T2 > &)
 Assign to this tensor without creating a temporary. More...
 
template<typename Scalar >
boostcopy::enable_if_c< ScalarTraits< Scalar >::value, TypeTensor & >::type operator= (const Scalar &libmesh_dbg_var(p))
 Assignment-from-scalar operator. More...
 
const T & operator() (const unsigned int i, const unsigned int j) const
 
T & operator() (const unsigned int i, const unsigned int j)
 
ConstTypeTensorColumn< T > slice (const unsigned int i) const
 
TypeTensorColumn< T > slice (const unsigned int i)
 
TypeVector< T > row (const unsigned int r) const
 
template<typename T2 >
TypeTensor< typename CompareTypes< T, T2 >::supertype > operator+ (const TypeTensor< T2 > &) const
 Add another tensor to this tensor. More...
 
template<typename T2 >
const TypeTensor< T > & operator+= (const TypeTensor< T2 > &)
 Add to this tensor. More...
 
template<typename T2 >
void add (const TypeTensor< T2 > &)
 Add to this tensor without creating a temporary. More...
 
template<typename T2 >
void add_scaled (const TypeTensor< T2 > &, const T)
 Add a scaled tensor to this tensor without creating a temporary. More...
 
template<typename T2 >
TypeTensor< typename CompareTypes< T, T2 >::supertype > operator- (const TypeTensor< T2 > &) const
 Subtract a tensor from this tensor. More...
 
template<typename T2 >
const TypeTensor< T > & operator-= (const TypeTensor< T2 > &)
 Subtract from this tensor. More...
 
template<typename T2 >
void subtract (const TypeTensor< T2 > &)
 Subtract from this tensor without creating a temporary. More...
 
template<typename T2 >
void subtract_scaled (const TypeTensor< T2 > &, const T)
 Subtract a scaled value from this tensor without creating a temporary. More...
 
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
 Multiply this tensor by a scalar value. More...
 
template<typename Scalar >
const TypeTensor< T > & operator*= (const Scalar)
 Multiply this tensor by a scalar value in place. More...
 
template<typename Scalar >
boostcopy::enable_if_c< ScalarTraits< Scalar >::value, TypeTensor< typename CompareTypes< T, Scalar >::supertype > >::type operator/ (const Scalar) const
 Divide each entry of this tensor by a scalar value. More...
 
const TypeTensor< T > & operator/= (const T)
 Divide each entry of this tensor by a scalar value. More...
 
template<typename T2 >
TypeTensor< T > operator* (const TypeTensor< T2 > &) const
 Multiply 2 tensors together, i.e. More...
 
template<typename T2 >
CompareTypes< T, T2 >::supertype contract (const TypeTensor< T2 > &) const
 Multiply 2 tensors together, i.e. More...
 
template<typename T2 >
TypeVector< typename CompareTypes< T, T2 >::supertype > operator* (const TypeVector< T2 > &) const
 Multiply this tensor by a vector, i.e. More...
 
TypeTensor< T > transpose () const
 
TypeTensor< T > inverse () const
 
void solve (const TypeVector< T > &b, TypeVector< T > &x) const
 Solve the 2x2 or 3x3 system of equations A*x = b for x by directly inverting A and multiplying it by b. More...
 
Real size () const
 
Real norm () const
 
Real size_sq () const
 
Real norm_sq () const
 
det () const
 
tr () const
 
void zero ()
 Set all entries of the tensor to 0. More...
 
bool operator== (const TypeTensor< T > &rhs) const
 
bool operator< (const TypeTensor< T > &rhs) const
 
bool operator> (const TypeTensor< T > &rhs) const
 
void print (std::ostream &os=libMesh::out) const
 Formatted print to a stream which defaults to libMesh::out. More...
 
void write_unformatted (std::ostream &out, const bool newline=true) const
 Unformatted print to the stream out. More...
 
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 ()
 Empty constructor. More...
 
 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)
 Constructor-from-T. More...
 
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)
 Constructor-from-Scalar. More...
 
template<typename T2 >
 TypeTensor (const TypeVector< T2 > &vx)
 Constructor. More...
 
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]
 The coordinates of the TypeTensor. More...
 

Friends

template<typename T2 >
class TypeTensor
 
std::ostream & operator<< (std::ostream &os, const TypeTensor< T > &t)
 Formatted print as above but supports the syntax: More...
 

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

Empty constructor.

Gives the tensor 0 in LIBMESH_DIM dimensions.

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

484 {
485  _coords[0] = 0;
486 
487 #if LIBMESH_DIM > 1
488  _coords[1] = 0;
489  _coords[2] = 0;
490  _coords[3] = 0;
491 #endif
492 
493 #if LIBMESH_DIM > 2
494  _coords[4] = 0;
495  _coords[5] = 0;
496  _coords[6] = 0;
497  _coords[7] = 0;
498  _coords[8] = 0;
499 #endif
500 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
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 
)
explicitprotected

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

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

515 {
516  _coords[0] = xx;
517 
518 #if LIBMESH_DIM == 2
519  _coords[1] = xy;
520  _coords[2] = yx;
521  _coords[3] = yy;
522 #endif
523 
524 #if LIBMESH_DIM == 3
525  _coords[1] = xy;
526  _coords[2] = xz;
527  _coords[3] = yx;
528  _coords[4] = yy;
529  _coords[5] = yz;
530  _coords[6] = zx;
531  _coords[7] = zy;
532  _coords[8] = zz;
533 #endif
534 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
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 
)
explicitprotected

Constructor-from-Scalar.

Definition at line 540 of file type_tensor.h.

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

551 {
552  _coords[0] = xx;
553 
554 #if LIBMESH_DIM == 2
555  _coords[1] = xy;
556  _coords[2] = yx;
557  _coords[3] = yy;
558 #endif
559 
560 #if LIBMESH_DIM == 3
561  _coords[1] = xy;
562  _coords[2] = xz;
563  _coords[3] = yx;
564  _coords[4] = yy;
565  _coords[5] = yz;
566  _coords[6] = zx;
567  _coords[7] = zy;
568  _coords[8] = zz;
569 #endif
570 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
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 587 of file type_tensor.h.

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

588 {
589  libmesh_assert_equal_to (LIBMESH_DIM, 1);
590  _coords[0] = vx(0);
591 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
template<typename T >
template<typename T2 >
libMesh::TypeTensor< T >::TypeTensor ( const TypeVector< T2 > &  vx,
const TypeVector< T2 > &  vy 
)
protected

Definition at line 595 of file type_tensor.h.

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

597 {
598  libmesh_assert_equal_to (LIBMESH_DIM, 2);
599  _coords[0] = vx(0);
600  _coords[1] = vx(1);
601  _coords[2] = vy(0);
602  _coords[3] = vy(1);
603 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
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 607 of file type_tensor.h.

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

610 {
611  libmesh_assert_equal_to (LIBMESH_DIM, 3);
612  _coords[0] = vx(0);
613  _coords[1] = vx(1);
614  _coords[2] = vx(2);
615  _coords[3] = vy(0);
616  _coords[4] = vy(1);
617  _coords[5] = vy(2);
618  _coords[6] = vz(0);
619  _coords[7] = vz(1);
620  _coords[8] = vz(2);
621 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
template<typename T >
template<typename T2 >
libMesh::TypeTensor< T >::TypeTensor ( const TypeTensor< T2 > &  p)

Copy-constructor.

Definition at line 577 of file type_tensor.h.

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

578 {
579  // copy the nodes from vector p to me
580  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
581  _coords[i] = p._coords[i];
582 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
template<typename T >
libMesh::TypeTensor< T >::~TypeTensor ( )

Destructor.

Definition at line 628 of file type_tensor.h.

629 {
630 }

Member Function Documentation

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

Add to this tensor without creating a temporary.

Definition at line 766 of file type_tensor.h.

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

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

767 {
768  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
769  _coords[i] += p._coords[i];
770 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
template<typename T >
template<typename T2 >
void libMesh::TypeTensor< T >::add_scaled ( const TypeTensor< T2 > &  p,
const T  factor 
)

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

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

778 {
779  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
780  _coords[i] += factor*p._coords[i];
781 
782 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
template<typename T >
template<typename T2 >
void libMesh::TypeTensor< T >::assign ( const TypeTensor< T2 > &  p)

Assign to this tensor without creating a temporary.

Definition at line 637 of file type_tensor.h.

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

638 {
639  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
640  _coords[i] = p._coords[i];
641 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
template<typename T >
template<typename T2 >
CompareTypes< T, T2 >::supertype libMesh::TypeTensor< T >::contract ( const TypeTensor< T2 > &  t) const

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.

sum Aij*Bij. The tensors may be of different types.

Definition at line 1175 of file type_tensor.h.

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

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

1176 {
1177  typename CompareTypes<T,T2>::supertype sum = 0.;
1178  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1179  sum += _coords[i]*t._coords[i];
1180  return sum;
1181 }
void sum(T &r, const Communicator &comm=Communicator_World)
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
template<typename T >
T libMesh::TypeTensor< T >::det ( ) const
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 1208 of file type_tensor.h.

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

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

1209 {
1210 #if LIBMESH_DIM == 1
1211  return _coords[0];
1212 #endif
1213 
1214 #if LIBMESH_DIM == 2
1215  return (_coords[0] * _coords[3]
1216  - _coords[1] * _coords[2]);
1217 #endif
1218 
1219 #if LIBMESH_DIM == 3
1220  return (_coords[0] * _coords[4] * _coords[8]
1221  + _coords[1] * _coords[5] * _coords[6]
1222  + _coords[2] * _coords[3] * _coords[7]
1223  - _coords[0] * _coords[5] * _coords[7]
1224  - _coords[1] * _coords[3] * _coords[8]
1225  - _coords[2] * _coords[4] * _coords[6]);
1226 #endif
1227 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
template<typename T >
TypeTensor< T > libMesh::TypeTensor< T >::inverse ( ) const
Returns
The inverse of this tensor as an independent object.

Definition at line 1022 of file type_tensor.h.

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

Referenced by NonlinearNeoHookeCurrentConfig::calculate_stress(), libMesh::TypeTensor< T >::operator=(), and TypeTensorTest::testInverse().

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

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

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

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

1271 {
1272  Real sum = 0.;
1273  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1274  sum += TensorTools::norm_sq(_coords[i]);
1275  return sum;
1276 }
void sum(T &r, const Communicator &comm=Communicator_World)
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
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
Returns
A const reference to the (i,j) entry of the tensor.

Definition at line 647 of file type_tensor.h.

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

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

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

Definition at line 666 of file type_tensor.h.

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

668 {
669 #if LIBMESH_DIM < 3
670 
671  if (i >= LIBMESH_DIM || j >= LIBMESH_DIM)
672  libmesh_error_msg("ERROR: You are assigning to a tensor component that is out of range for the compiled LIBMESH_DIM!");
673 
674 #endif
675 
676  libmesh_assert_less (i, LIBMESH_DIM);
677  libmesh_assert_less (j, LIBMESH_DIM);
678 
679  return _coords[i*LIBMESH_DIM+j];
680 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
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

Multiply this tensor by a scalar value.

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

Definition at line 893 of file type_tensor.h.

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

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

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

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

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

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

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

Multiply this tensor by a scalar value in place.

Returns
A reference to *this.

Definition at line 940 of file type_tensor.h.

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

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

941 {
942  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
943  _coords[i] *= factor;
944 
945  return *this;
946 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
template<typename T >
template<typename T2 >
TypeTensor< typename CompareTypes< T, T2 >::supertype > libMesh::TypeTensor< T >::operator+ ( const TypeTensor< T2 > &  p) const

Add another tensor to this tensor.

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

Definition at line 720 of file type_tensor.h.

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

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

721 {
722 
723 #if LIBMESH_DIM == 1
724  return TypeTensor(_coords[0] + p._coords[0]);
725 #endif
726 
727 #if LIBMESH_DIM == 2
728  return TypeTensor(_coords[0] + p._coords[0],
729  _coords[1] + p._coords[1],
730  0.,
731  _coords[2] + p._coords[2],
732  _coords[3] + p._coords[3]);
733 #endif
734 
735 #if LIBMESH_DIM == 3
736  return TypeTensor(_coords[0] + p._coords[0],
737  _coords[1] + p._coords[1],
738  _coords[2] + p._coords[2],
739  _coords[3] + p._coords[3],
740  _coords[4] + p._coords[4],
741  _coords[5] + p._coords[5],
742  _coords[6] + p._coords[6],
743  _coords[7] + p._coords[7],
744  _coords[8] + p._coords[8]);
745 #endif
746 
747 }
TypeTensor()
Empty constructor.
Definition: type_tensor.h:483
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
template<typename T >
template<typename T2 >
const TypeTensor< T > & libMesh::TypeTensor< T >::operator+= ( const TypeTensor< T2 > &  p)

Add to this tensor.

Returns
A reference to *this.

Definition at line 754 of file type_tensor.h.

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

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

755 {
756  this->add (p);
757 
758  return *this;
759 }
void add(const TypeTensor< T2 > &)
Add to this tensor without creating a temporary.
Definition: type_tensor.h:766
template<typename T >
template<typename T2 >
TypeTensor< typename CompareTypes< T, T2 >::supertype > libMesh::TypeTensor< T >::operator- ( const TypeTensor< T2 > &  p) const

Subtract a tensor from this tensor.

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

Definition at line 790 of file type_tensor.h.

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

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

Definition at line 857 of file type_tensor.h.

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

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

858 {
859 
860 #if LIBMESH_DIM == 1
861  return TypeTensor(-_coords[0]);
862 #endif
863 
864 #if LIBMESH_DIM == 2
865  return TypeTensor(-_coords[0],
866  -_coords[1],
867  -_coords[2],
868  -_coords[3]);
869 #endif
870 
871 #if LIBMESH_DIM == 3
872  return TypeTensor(-_coords[0],
873  -_coords[1],
874  -_coords[2],
875  -_coords[3],
876  -_coords[4],
877  -_coords[5],
878  -_coords[6],
879  -_coords[7],
880  -_coords[8]);
881 #endif
882 
883 }
TypeTensor()
Empty constructor.
Definition: type_tensor.h:483
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
template<typename T >
template<typename T2 >
const TypeTensor< T > & libMesh::TypeTensor< T >::operator-= ( const TypeTensor< T2 > &  p)

Subtract from this tensor.

Returns
A reference to *this.

Definition at line 824 of file type_tensor.h.

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

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

825 {
826  this->subtract (p);
827 
828  return *this;
829 }
void subtract(const TypeTensor< T2 > &)
Subtract from this tensor without creating a temporary.
Definition: type_tensor.h:836
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

Divide each entry of this tensor by a scalar value.

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

Definition at line 957 of file type_tensor.h.

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

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

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

Divide each entry of this tensor by a scalar value.

Returns
A reference to *this.

Definition at line 1122 of file type_tensor.h.

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

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

1123 {
1124  libmesh_assert_not_equal_to (factor, static_cast<T>(0.));
1125 
1126  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1127  _coords[i] /= factor;
1128 
1129  return *this;
1130 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
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)
template<typename T >
bool libMesh::TypeTensor< T >::operator== ( const TypeTensor< T > &  rhs) const
Returns
true if two tensors are equal, false otherwise.

Definition at line 1282 of file type_tensor.h.

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

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

1283 {
1284 #if LIBMESH_DIM == 1
1285  return (std::abs(_coords[0] - rhs._coords[0])
1286  < TOLERANCE);
1287 #endif
1288 
1289 #if LIBMESH_DIM == 2
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  < 4.*TOLERANCE);
1295 #endif
1296 
1297 #if LIBMESH_DIM == 3
1298  return ((std::abs(_coords[0] - rhs._coords[0]) +
1299  std::abs(_coords[1] - rhs._coords[1]) +
1300  std::abs(_coords[2] - rhs._coords[2]) +
1301  std::abs(_coords[3] - rhs._coords[3]) +
1302  std::abs(_coords[4] - rhs._coords[4]) +
1303  std::abs(_coords[5] - rhs._coords[5]) +
1304  std::abs(_coords[6] - rhs._coords[6]) +
1305  std::abs(_coords[7] - rhs._coords[7]) +
1306  std::abs(_coords[8] - rhs._coords[8]))
1307  < 9.*TOLERANCE);
1308 #endif
1309 
1310 }
double abs(double a)
static const Real TOLERANCE
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
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
Returns
A copy of one row of the tensor as a TypeVector.

Definition at line 706 of file type_tensor.h.

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

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

707 {
708  TypeVector<T> return_vector;
709 
710  for (unsigned int j=0; j<LIBMESH_DIM; j++)
711  return_vector._coords[j] = _coords[r*LIBMESH_DIM + j];
712 
713  return return_vector;
714 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
template<typename T >
Real libMesh::TypeTensor< T >::size ( ) const
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 1188 of file type_tensor.h.

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

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

1189 {
1190  libmesh_deprecated();
1191  return this->norm();
1192 }
Real norm() const
Definition: type_tensor.h:1199
template<typename T >
Real libMesh::TypeTensor< T >::size_sq ( ) const
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 1259 of file type_tensor.h.

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

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

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

Definition at line 686 of file type_tensor.h.

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

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

Definition at line 696 of file type_tensor.h.

697 {
698  libmesh_assert_less (i, LIBMESH_DIM);
699  return TypeTensorColumn<T>(*this, i);
700 }
template<typename T >
void libMesh::TypeTensor< T >::solve ( const TypeVector< T > &  b,
TypeVector< T > &  x 
) const

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

Returns
The solution in the x vector.

Definition at line 1071 of file type_tensor.h.

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

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

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

Subtract from this tensor without creating a temporary.

Definition at line 836 of file type_tensor.h.

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

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

837 {
838  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
839  _coords[i] -= p._coords[i];
840 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
template<typename T >
template<typename T2 >
void libMesh::TypeTensor< T >::subtract_scaled ( const TypeTensor< T2 > &  p,
const T  factor 
)

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

Definition at line 847 of file type_tensor.h.

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

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

848 {
849  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
850  _coords[i] -= factor*p._coords[i];
851 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
template<typename T >
T libMesh::TypeTensor< T >::tr ( ) const
Returns
The trace of the tensor.

Definition at line 1231 of file type_tensor.h.

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

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

1232 {
1233 #if LIBMESH_DIM == 1
1234  return _coords[0];
1235 #endif
1236 
1237 #if LIBMESH_DIM == 2
1238  return _coords[0] + _coords[3];
1239 #endif
1240 
1241 #if LIBMESH_DIM == 3
1242  return _coords[0] + _coords[4] + _coords[8];
1243 #endif
1244 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
template<typename T >
TypeTensor< T > libMesh::TypeTensor< T >::transpose ( ) const
Returns
The transpose of this tensor (with complex numbers not conjugated).

Definition at line 992 of file type_tensor.h.

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

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

993 {
994 #if LIBMESH_DIM == 1
995  return TypeTensor(_coords[0]);
996 #endif
997 
998 #if LIBMESH_DIM == 2
999  return TypeTensor(_coords[0],
1000  _coords[2],
1001  _coords[1],
1002  _coords[3]);
1003 #endif
1004 
1005 #if LIBMESH_DIM == 3
1006  return TypeTensor(_coords[0],
1007  _coords[3],
1008  _coords[6],
1009  _coords[1],
1010  _coords[4],
1011  _coords[7],
1012  _coords[2],
1013  _coords[5],
1014  _coords[8]);
1015 #endif
1016 }
TypeTensor()
Empty constructor.
Definition: type_tensor.h:483
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417
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 ( )

Set all entries of the tensor to 0.

Definition at line 1248 of file type_tensor.h.

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

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

1249 {
1250  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1251  _coords[i] = 0.;
1252 }
T _coords[LIBMESH_DIM *LIBMESH_DIM]
The coordinates of the TypeTensor.
Definition: type_tensor.h:417

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

401  {
402  t.print(os);
403  return os;
404  }
template<typename T>
template<typename T2 >
friend class TypeTensor
friend

Definition at line 52 of file type_tensor.h.

Member Data Documentation

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

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