libMesh::TensorValue< T > Class Template Reference

#include <exact_error_estimator.h>

Inheritance diagram for libMesh::TensorValue< T >:

Public Member Functions

 TensorValue ()
 
 TensorValue (const T xx, const T xy=0, const T xz=0, const T yx=0, const T yy=0, const T yz=0, const T zx=0, const T zy=0, const T zz=0)
 
template<typename Scalar >
 TensorValue (const Scalar xx, const Scalar xy=0, const Scalar xz=0, const Scalar yx=0, const Scalar yy=0, const Scalar yz=0, const Scalar zx=0, const Scalar zy=0, typename boostcopy::enable_if_c< ScalarTraits< Scalar >::value, const Scalar >::type zz=0)
 
template<typename T2 >
 TensorValue (const TypeVector< T2 > &vx)
 
template<typename T2 >
 TensorValue (const TypeVector< T2 > &vx, const TypeVector< T2 > &vy)
 
template<typename T2 >
 TensorValue (const TypeVector< T2 > &vx, const TypeVector< T2 > &vy, const TypeVector< T2 > &vz)
 
template<typename T2 >
 TensorValue (const TensorValue< T2 > &p)
 
template<typename T2 >
 TensorValue (const TypeTensor< T2 > &p)
 
 TensorValue (const TypeTensor< Real > &p_re, const TypeTensor< Real > &p_im)
 
template<typename Scalar >
boostcopy::enable_if_c
< ScalarTraits< Scalar >
::value, TensorValue & >::type 
operator= (const Scalar &p)
 
template<typename T2 >
void assign (const TypeTensor< T2 > &)
 
const T & operator() (const unsigned int i, const unsigned int j) const
 
T & operator() (const unsigned int i, const unsigned int j)
 
ConstTypeTensorColumn< T > slice (const unsigned int i) const
 
TypeTensorColumn< T > slice (const unsigned int i)
 
TypeVector< T > row (const unsigned int r)
 
template<typename T2 >
TypeTensor< typename
CompareTypes< T, T2 >
::supertype > 
operator+ (const TypeTensor< T2 > &) const
 
template<typename T2 >
const TypeTensor< T > & operator+= (const TypeTensor< T2 > &)
 
template<typename T2 >
void add (const TypeTensor< T2 > &)
 
template<typename T2 >
void add_scaled (const TypeTensor< T2 > &, const T)
 
template<typename T2 >
TypeTensor< typename
CompareTypes< T, T2 >
::supertype > 
operator- (const TypeTensor< T2 > &) const
 
TypeTensor< T > operator- () const
 
template<typename T2 >
const TypeTensor< T > & operator-= (const TypeTensor< T2 > &)
 
template<typename T2 >
void subtract (const TypeTensor< T2 > &)
 
template<typename T2 >
void subtract_scaled (const TypeTensor< T2 > &, const T)
 
template<typename Scalar >
boostcopy::enable_if_c
< ScalarTraits< Scalar >
::value, TypeTensor< typename
CompareTypes< T, Scalar >
::supertype > >::type 
operator* (const Scalar) const
 
template<typename T2 >
TypeTensor< T > operator* (const TypeTensor< T2 > &) const
 
template<typename T2 >
TypeVector< typename
CompareTypes< T, T2 >
::supertype > 
operator* (const TypeVector< T2 > &) const
 
template<typename Scalar >
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 >
CompareTypes< T, T2 >::supertype contract (const TypeTensor< T2 > &) const
 
TypeTensor< T > transpose () const
 
Real size () const
 
Real size_sq () const
 
det () const
 
tr () const
 
void zero ()
 
bool operator== (const TypeTensor< T > &rhs) const
 
bool operator< (const TypeTensor< T > &rhs) const
 
template<>
bool operator< (const TypeTensor< Real > &rhs) const
 
template<>
bool operator< (const TypeTensor< Complex > &rhs) const
 
bool operator> (const TypeTensor< T > &rhs) const
 
template<>
bool operator> (const TypeTensor< Real > &rhs) const
 
template<>
bool operator> (const TypeTensor< Complex > &rhs) const
 
void print (std::ostream &os=libMesh::out) const
 
void write_unformatted (std::ostream &out, const bool newline=true) const
 

Public Attributes

_coords [LIBMESH_DIM *LIBMESH_DIM]
 

Detailed Description

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

This class defines a tensor in LIBMESH_DIM dimensional Real or Complex space. The typedef RealTensorValue always defines a real-valued tensor, and NumberTensorValue defines a real or complex-valued tensor depending on how the library was configured.

Author
Roy H. Stogner 2004

Definition at line 46 of file exact_error_estimator.h.

Constructor & Destructor Documentation

template<typename T >
libMesh::TensorValue< T >::TensorValue ( )
inline

Empty constructor. Gives the tensor 0 in LIBMESH_DIM dimensional T space.

Definition at line 157 of file tensor_value.h.

157  :
158  TypeTensor<T> ()
159 {
160 }
template<typename T >
libMesh::TensorValue< T >::TensorValue ( const T  xx,
const T  xy = 0,
const T  xz = 0,
const T  yx = 0,
const T  yy = 0,
const T  yz = 0,
const T  zx = 0,
const T  zy = 0,
const T  zz = 0 
)
inlineexplicit

Constructor-from-T. By default sets higher dimensional entries to 0.

Definition at line 167 of file tensor_value.h.

175  :
176  TypeTensor<T> (xx,xy,xz,yx,yy,yz,zx,zy,zz)
177 {
178 }
template<typename T >
template<typename Scalar >
libMesh::TensorValue< T >::TensorValue ( const Scalar  xx,
const Scalar  xy = 0,
const Scalar  xz = 0,
const Scalar  yx = 0,
const Scalar  yy = 0,
const Scalar  yz = 0,
const Scalar  zx = 0,
const Scalar  zy = 0,
typename boostcopy::enable_if_c< ScalarTraits< Scalar >::value, const Scalar >::type  zz = 0 
)
inlineexplicit

Constructor-from-scalars. By default sets higher dimensional entries to 0.

Definition at line 185 of file tensor_value.h.

195  :
196  TypeTensor<T> (xx,xy,xz,yx,yy,yz,zx,zy,zz)
197 {
198 }
template<typename T >
template<typename T2 >
libMesh::TensorValue< T >::TensorValue ( const TypeVector< T2 > &  vx)
inline

Constructor. Takes 1 row vector for LIBMESH_DIM=1

Definition at line 215 of file tensor_value.h.

215  :
216  TypeTensor<T> (vx)
217 {
218 }
template<typename T >
template<typename T2 >
libMesh::TensorValue< T >::TensorValue ( const TypeVector< T2 > &  vx,
const TypeVector< T2 > &  vy 
)
inline

Constructor. Takes 2 row vectors for LIBMESH_DIM=2

Definition at line 225 of file tensor_value.h.

226  :
227  TypeTensor<T> (vx, vy)
228 {
229 }
template<typename T >
template<typename T2 >
libMesh::TensorValue< T >::TensorValue ( const TypeVector< T2 > &  vx,
const TypeVector< T2 > &  vy,
const TypeVector< T2 > &  vz 
)
inline

Constructor. Takes 3 row vectors for LIBMESH_DIM=3

Definition at line 236 of file tensor_value.h.

238  :
239  TypeTensor<T> (vx, vy, vz)
240 {
241 }
template<typename T >
template<typename T2 >
libMesh::TensorValue< T >::TensorValue ( const TensorValue< T2 > &  p)
inline

Copy-constructor.

Definition at line 205 of file tensor_value.h.

205  :
206  TypeTensor<T> (p)
207 {
208 }
template<typename T >
template<typename T2 >
libMesh::TensorValue< T >::TensorValue ( const TypeTensor< T2 > &  p)
inline

Copy-constructor.

Definition at line 248 of file tensor_value.h.

248  :
249  TypeTensor<T> (p)
250 {
251 }
template<typename T >
libMesh::TensorValue< T >::TensorValue ( const TypeTensor< Real > &  p_re,
const TypeTensor< Real > &  p_im 
)
inline

Constructor that takes two TypeTensor<Real> representing the real and imaginary part as arguments.

Definition at line 257 of file tensor_value.h.

258  :
259  TypeTensor<T> (Complex (p_re(0,0), p_im(0,0)),
260  Complex (p_re(0,1), p_im(0,1)),
261  Complex (p_re(0,2), p_im(0,2)),
262  Complex (p_re(1,0), p_im(1,0)),
263  Complex (p_re(1,1), p_im(1,1)),
264  Complex (p_re(1,2), p_im(1,2)),
265  Complex (p_re(2,0), p_im(2,0)),
266  Complex (p_re(2,1), p_im(2,1)),
267  Complex (p_re(2,2), p_im(2,2)))
268 {
269 }

Member Function Documentation

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

Add to this tensor without creating a temporary.

Definition at line 706 of file type_tensor.h.

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

707 {
708  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
709  _coords[i] += p._coords[i];
710 }
template<typename T >
template<typename T2 >
void libMesh::TypeTensor< T >::add_scaled ( const TypeTensor< T2 > &  p,
const T  factor 
)
inlineinherited

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

Definition at line 717 of file type_tensor.h.

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

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

718 {
719  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
720  _coords[i] += factor*p._coords[i];
721 
722 }
template<typename T >
template<typename T2 >
void libMesh::TypeTensor< T >::assign ( const TypeTensor< T2 > &  p)
inlineinherited

Assign to a tensor without creating a temporary.

Definition at line 572 of file type_tensor.h.

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

573 {
574  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
575  _coords[i] = p._coords[i];
576 }
template<typename T >
template<typename T2 >
CompareTypes< T, T2 >::supertype libMesh::TypeTensor< T >::contract ( const TypeTensor< T2 > &  t) const
inlineinherited

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

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

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

1016 {
1017  typename CompareTypes<T,T2>::supertype sum = 0.;
1018  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1019  sum += _coords[i]*t._coords[i];
1020  return sum;
1021 }
template<typename T >
T libMesh::TypeTensor< T >::det ( ) const
inlineinherited

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

Definition at line 1036 of file type_tensor.h.

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

1037 {
1038 #if LIBMESH_DIM == 1
1039  return _coords[0];
1040 #endif
1041 
1042 #if LIBMESH_DIM == 2
1043  return (_coords[0] * _coords[3]
1044  - _coords[1] * _coords[2]);
1045 #endif
1046 
1047 #if LIBMESH_DIM == 3
1048  return (_coords[0] * _coords[4] * _coords[8]
1049  + _coords[1] * _coords[5] * _coords[6]
1050  + _coords[2] * _coords[3] * _coords[7]
1051  - _coords[0] * _coords[5] * _coords[7]
1052  - _coords[1] * _coords[3] * _coords[8]
1053  - _coords[2] * _coords[4] * _coords[6]);
1054 #endif
1055 }
template<typename T >
const T & libMesh::TypeTensor< T >::operator() ( const unsigned int  i,
const unsigned int  j 
) const
inlineinherited

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

Definition at line 582 of file type_tensor.h.

584 {
585  libmesh_assert_less (i, 3);
586  libmesh_assert_less (j, 3);
587 
588 #if LIBMESH_DIM < 3
589  const static T my_zero = 0;
590  if (i >= LIBMESH_DIM || j >= LIBMESH_DIM)
591  return my_zero;
592 #endif
593 
594  return _coords[i*LIBMESH_DIM+j];
595 }
template<typename T >
T & libMesh::TypeTensor< T >::operator() ( const unsigned int  i,
const unsigned int  j 
)
inlineinherited

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

Definition at line 601 of file type_tensor.h.

603 {
604 #if LIBMESH_DIM < 3
605 
606  if (i >= LIBMESH_DIM || j >= LIBMESH_DIM)
607  {
608 // libMesh::err << "ERROR: You are assigning to a tensor component" << std::endl
609 // << "that is out of range for the compiled LIBMESH_DIM!" << std::endl
610 // << " LIBMESH_DIM=" << LIBMESH_DIM << " , i=" << i << " , j=" << j << std::endl;
611  libmesh_error();
612  }
613 
614 #endif
615 
616  libmesh_assert_less (i, LIBMESH_DIM);
617  libmesh_assert_less (j, LIBMESH_DIM);
618 
619  return _coords[i*LIBMESH_DIM+j];
620 }
template<typename T >
template<typename Scalar >
boostcopy::enable_if_c< ScalarTraits< Scalar >::value, TypeTensor< typename CompareTypes< T, Scalar >::supertype > >::type libMesh::TypeTensor< T >::operator* ( const Scalar  factor) const
inlineinherited

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

Definition at line 833 of file type_tensor.h.

834 {
835  typedef typename CompareTypes<T, Scalar>::supertype TS;
836 
837 
838 #if LIBMESH_DIM == 1
839  return TypeTensor<TS>(_coords[0]*factor);
840 #endif
841 
842 #if LIBMESH_DIM == 2
843  return TypeTensor<TS>(_coords[0]*factor,
844  _coords[1]*factor,
845  _coords[2]*factor,
846  _coords[3]*factor);
847 #endif
848 
849 #if LIBMESH_DIM == 3
850  return TypeTensor<TS>(_coords[0]*factor,
851  _coords[1]*factor,
852  _coords[2]*factor,
853  _coords[3]*factor,
854  _coords[4]*factor,
855  _coords[5]*factor,
856  _coords[6]*factor,
857  _coords[7]*factor,
858  _coords[8]*factor);
859 #endif
860 }
template<typename T >
template<typename T2 >
TypeTensor< T > libMesh::TypeTensor< T >::operator* ( const TypeTensor< T2 > &  p) const
inlineinherited

Multiply 2 tensors together, i.e. matrix product. The tensors may be of different types.

Definition at line 994 of file type_tensor.h.

995 {
996  TypeTensor<T> returnval;
997  for (unsigned int i=0; i<LIBMESH_DIM; i++)
998  for (unsigned int j=0; j<LIBMESH_DIM; j++)
999  for (unsigned int k=0; k<LIBMESH_DIM; k++)
1000  returnval(i,j) += (*this)(i,k)*p(k,j);
1001 
1002  return returnval;
1003 }
template<typename T >
template<typename T2 >
TypeVector< typename CompareTypes< T, T2 >::supertype > libMesh::TypeTensor< T >::operator* ( const TypeVector< T2 > &  p) const
inlineinherited

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

Definition at line 979 of file type_tensor.h.

980 {
981  TypeVector<typename CompareTypes<T,T2>::supertype> returnval;
982  for (unsigned int i=0; i<LIBMESH_DIM; i++)
983  for (unsigned int j=0; j<LIBMESH_DIM; j++)
984  returnval(i) += (*this)(i,j)*p(j);
985 
986  return returnval;
987 }
template<typename T >
template<typename Scalar >
const TypeTensor< T > & libMesh::TypeTensor< T >::operator*= ( const Scalar  factor)
inlineinherited

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

Definition at line 880 of file type_tensor.h.

881 {
882  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
883  _coords[i] *= factor;
884 
885  return *this;
886 }
template<typename T >
template<typename T2 >
TypeTensor< typename CompareTypes< T, T2 >::supertype > libMesh::TypeTensor< T >::operator+ ( const TypeTensor< T2 > &  p) const
inlineinherited

Add two tensors.

Definition at line 660 of file type_tensor.h.

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

661 {
662 
663 #if LIBMESH_DIM == 1
664  return TypeTensor(_coords[0] + p._coords[0]);
665 #endif
666 
667 #if LIBMESH_DIM == 2
668  return TypeTensor(_coords[0] + p._coords[0],
669  _coords[1] + p._coords[1],
670  0.,
671  _coords[2] + p._coords[2],
672  _coords[3] + p._coords[3]);
673 #endif
674 
675 #if LIBMESH_DIM == 3
676  return TypeTensor(_coords[0] + p._coords[0],
677  _coords[1] + p._coords[1],
678  _coords[2] + p._coords[2],
679  _coords[3] + p._coords[3],
680  _coords[4] + p._coords[4],
681  _coords[5] + p._coords[5],
682  _coords[6] + p._coords[6],
683  _coords[7] + p._coords[7],
684  _coords[8] + p._coords[8]);
685 #endif
686 
687 }
template<typename T >
template<typename T2 >
const TypeTensor< T > & libMesh::TypeTensor< T >::operator+= ( const TypeTensor< T2 > &  p)
inlineinherited

Add to this tensor.

Definition at line 694 of file type_tensor.h.

695 {
696  this->add (p);
697 
698  return *this;
699 }
template<typename T >
template<typename T2 >
TypeTensor< typename CompareTypes< T, T2 >::supertype > libMesh::TypeTensor< T >::operator- ( const TypeTensor< T2 > &  p) const
inlineinherited

Subtract two tensors.

Definition at line 730 of file type_tensor.h.

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

731 {
732 
733 #if LIBMESH_DIM == 1
734  return TypeTensor(_coords[0] - p._coords[0]);
735 #endif
736 
737 #if LIBMESH_DIM == 2
738  return TypeTensor(_coords[0] - p._coords[0],
739  _coords[1] - p._coords[1],
740  0.,
741  _coords[2] - p._coords[2],
742  _coords[3] - p._coords[3]);
743 #endif
744 
745 #if LIBMESH_DIM == 3
746  return TypeTensor(_coords[0] - p._coords[0],
747  _coords[1] - p._coords[1],
748  _coords[2] - p._coords[2],
749  _coords[3] - p._coords[3],
750  _coords[4] - p._coords[4],
751  _coords[5] - p._coords[5],
752  _coords[6] - p._coords[6],
753  _coords[7] - p._coords[7],
754  _coords[8] - p._coords[8]);
755 #endif
756 
757 }
template<typename T >
TypeTensor< T > libMesh::TypeTensor< T >::operator- ( ) const
inlineinherited

Return the opposite of a tensor

Definition at line 797 of file type_tensor.h.

798 {
799 
800 #if LIBMESH_DIM == 1
801  return TypeTensor(-_coords[0]);
802 #endif
803 
804 #if LIBMESH_DIM == 2
805  return TypeTensor(-_coords[0],
806  -_coords[1],
807  -_coords[2],
808  -_coords[3]);
809 #endif
810 
811 #if LIBMESH_DIM == 3
812  return TypeTensor(-_coords[0],
813  -_coords[1],
814  -_coords[2],
815  -_coords[3],
816  -_coords[4],
817  -_coords[5],
818  -_coords[6],
819  -_coords[7],
820  -_coords[8]);
821 #endif
822 
823 }
template<typename T >
template<typename T2 >
const TypeTensor< T > & libMesh::TypeTensor< T >::operator-= ( const TypeTensor< T2 > &  p)
inlineinherited

Subtract from this tensor.

Definition at line 764 of file type_tensor.h.

765 {
766  this->subtract (p);
767 
768  return *this;
769 }
template<typename T >
template<typename Scalar >
boostcopy::enable_if_c< ScalarTraits< Scalar >::value, TypeTensor< typename CompareTypes< T, Scalar >::supertype > >::type libMesh::TypeTensor< T >::operator/ ( const Scalar  factor) const
inlineinherited

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

Definition at line 897 of file type_tensor.h.

898 {
899  libmesh_assert_not_equal_to (factor, static_cast<T>(0.));
900 
901  typedef typename CompareTypes<T, Scalar>::supertype TS;
902 
903 #if LIBMESH_DIM == 1
904  return TypeTensor<TS>(_coords[0]/factor);
905 #endif
906 
907 #if LIBMESH_DIM == 2
908  return TypeTensor<TS>(_coords[0]/factor,
909  _coords[1]/factor,
910  _coords[2]/factor,
911  _coords[3]/factor);
912 #endif
913 
914 #if LIBMESH_DIM == 3
915  return TypeTensor<TS>(_coords[0]/factor,
916  _coords[1]/factor,
917  _coords[2]/factor,
918  _coords[3]/factor,
919  _coords[4]/factor,
920  _coords[5]/factor,
921  _coords[6]/factor,
922  _coords[7]/factor,
923  _coords[8]/factor);
924 #endif
925 
926 }
template<typename T >
const TypeTensor< T > & libMesh::TypeTensor< T >::operator/= ( const T  factor)
inlineinherited

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

Definition at line 962 of file type_tensor.h.

963 {
964  libmesh_assert_not_equal_to (factor, static_cast<T>(0.));
965 
966  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
967  _coords[i] /= factor;
968 
969  return *this;
970 }
template<>
bool libMesh::TypeTensor< Real >::operator< ( const TypeTensor< Real > &  rhs) const
inherited

Definition at line 113 of file type_tensor.C.

114 {
115  for (unsigned int i=0; i<LIBMESH_DIM; i++)
116  for (unsigned int j=0; j<LIBMESH_DIM; j++)
117  {
118  if ((*this)(i,j) < rhs(i,j))
119  return true;
120  if ((*this)(i,j) > rhs(i,j))
121  return false;
122  }
123  return false;
124 }
template<>
bool libMesh::TypeTensor< Complex >::operator< ( const TypeTensor< Complex > &  rhs) const
inherited

Definition at line 146 of file type_tensor.C.

147 {
148  for (unsigned int i=0; i<LIBMESH_DIM; i++)
149  for (unsigned int j=0; j<LIBMESH_DIM; j++)
150  {
151  if ((*this)(i,j).real() < rhs(i,j).real())
152  return true;
153  if ((*this)(i,j).real() > rhs(i,j).real())
154  return false;
155  if ((*this)(i,j).imag() < rhs(i,j).imag())
156  return true;
157  if ((*this)(i,j).imag() > rhs(i,j).imag())
158  return false;
159  }
160  return false;
161 }
template<typename T>
bool libMesh::TypeTensor< T >::operator< ( const TypeTensor< T > &  rhs) const
inherited
Returns
true if this tensor is "less" than another. Useful for sorting.
template<typename T>
template<typename Scalar >
boostcopy::enable_if_c< ScalarTraits<Scalar>::value, TensorValue&>::type libMesh::TensorValue< T >::operator= ( const Scalar &  p)
inline

Assignment-from-scalar operator. Used only to zero out tensors.

Definition at line 131 of file tensor_value.h.

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

132  { libmesh_assert_equal_to (p, Scalar(0)); this->zero(); return *this; }
template<typename T >
bool libMesh::TypeTensor< T >::operator== ( const TypeTensor< T > &  rhs) const
inlineinherited
Returns
true if two tensors are equal valued.

Definition at line 1098 of file type_tensor.h.

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

1099 {
1100 #if LIBMESH_DIM == 1
1101  return (std::abs(_coords[0] - rhs._coords[0])
1102  < TOLERANCE);
1103 #endif
1104 
1105 #if LIBMESH_DIM == 2
1106  return ((std::abs(_coords[0] - rhs._coords[0]) +
1107  std::abs(_coords[1] - rhs._coords[1]) +
1108  std::abs(_coords[2] - rhs._coords[2]) +
1109  std::abs(_coords[3] - rhs._coords[3]))
1110  < 4.*TOLERANCE);
1111 #endif
1112 
1113 #if LIBMESH_DIM == 3
1114  return ((std::abs(_coords[0] - rhs._coords[0]) +
1115  std::abs(_coords[1] - rhs._coords[1]) +
1116  std::abs(_coords[2] - rhs._coords[2]) +
1117  std::abs(_coords[3] - rhs._coords[3]) +
1118  std::abs(_coords[4] - rhs._coords[4]) +
1119  std::abs(_coords[5] - rhs._coords[5]) +
1120  std::abs(_coords[6] - rhs._coords[6]) +
1121  std::abs(_coords[7] - rhs._coords[7]) +
1122  std::abs(_coords[8] - rhs._coords[8]))
1123  < 9.*TOLERANCE);
1124 #endif
1125 
1126 }
template<>
bool libMesh::TypeTensor< Real >::operator> ( const TypeTensor< Real > &  rhs) const
inherited

Definition at line 129 of file type_tensor.C.

130 {
131  for (unsigned int i=0; i<LIBMESH_DIM; i++)
132  for (unsigned int j=0; j<LIBMESH_DIM; j++)
133  {
134  if ((*this)(i,j) > rhs(i,j))
135  return true;
136  if ((*this)(i,j) < rhs(i,j))
137  return false;
138  }
139  return false;
140 }
template<>
bool libMesh::TypeTensor< Complex >::operator> ( const TypeTensor< Complex > &  rhs) const
inherited

Definition at line 166 of file type_tensor.C.

167 {
168  for (unsigned int i=0; i<LIBMESH_DIM; i++)
169  for (unsigned int j=0; j<LIBMESH_DIM; j++)
170  {
171  if ((*this)(i,j).real() > rhs(i,j).real())
172  return true;
173  if ((*this)(i,j).real() < rhs(i,j).real())
174  return false;
175  if ((*this)(i,j).imag() > rhs(i,j).imag())
176  return true;
177  if ((*this)(i,j).imag() < rhs(i,j).imag())
178  return false;
179  }
180  return false;
181 }
template<typename T>
bool libMesh::TypeTensor< T >::operator> ( const TypeTensor< T > &  rhs) const
inherited
Returns
true if this tensor is "greater" than another.
template<typename T >
void libMesh::TypeTensor< T >::print ( std::ostream &  os = libMesh::out) const
inherited

Formatted print, by default to libMesh::out.

Definition at line 39 of file type_tensor.C.

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

Return one row of the tensor as a TypeVector.

Definition at line 646 of file type_tensor.h.

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

647 {
648  TypeVector<T> return_vector;
649 
650  for(unsigned int j=0; j<LIBMESH_DIM; j++)
651  return_vector._coords[j] = _coords[r*LIBMESH_DIM + j];
652 
653  return return_vector;
654 }
template<typename T >
Real libMesh::TypeTensor< T >::size ( ) const
inlineinherited

Returns the Frobenius norm of the tensor, i.e. the square-root of the sum of the elements squared.

Definition at line 1027 of file type_tensor.h.

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

1028 {
1029  return std::sqrt(this->size_sq());
1030 }
template<typename T >
Real libMesh::TypeTensor< T >::size_sq ( ) const
inlineinherited

Returns the Frobenius norm of the tensor squared, i.e. sum of the element magnitudes squared.

Definition at line 1086 of file type_tensor.h.

References 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(), and libMesh::HPCoarsenTest::select_refinement().

1087 {
1088  Real sum = 0.;
1089  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1090  sum += TensorTools::norm_sq(_coords[i]);
1091  return sum;
1092 }
template<typename T >
ConstTypeTensorColumn< T > libMesh::TypeTensor< T >::slice ( const unsigned int  i) const
inlineinherited

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

Definition at line 626 of file type_tensor.h.

627 {
628  libmesh_assert_less (i, LIBMESH_DIM);
629  return ConstTypeTensorColumn<T>(*this, i);
630 }
template<typename T >
TypeTensorColumn< T > libMesh::TypeTensor< T >::slice ( const unsigned int  i)
inlineinherited

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

Definition at line 636 of file type_tensor.h.

637 {
638  libmesh_assert_less (i, LIBMESH_DIM);
639  return TypeTensorColumn<T>(*this, i);
640 }
template<typename T >
template<typename T2 >
void libMesh::TypeTensor< T >::subtract ( const TypeTensor< T2 > &  p)
inlineinherited

Subtract from this tensor without creating a temporary.

Definition at line 776 of file type_tensor.h.

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

777 {
778  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
779  _coords[i] -= p._coords[i];
780 }
template<typename T >
template<typename T2 >
void libMesh::TypeTensor< T >::subtract_scaled ( const TypeTensor< T2 > &  p,
const T  factor 
)
inlineinherited

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

Definition at line 787 of file type_tensor.h.

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

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

788 {
789  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
790  _coords[i] -= factor*p._coords[i];
791 }
template<typename T >
T libMesh::TypeTensor< T >::tr ( ) const
inlineinherited

Returns the trace of the tensor.

Definition at line 1059 of file type_tensor.h.

1060 {
1061 #if LIBMESH_DIM == 1
1062  return _coords[0];
1063 #endif
1064 
1065 #if LIBMESH_DIM == 2
1066  return _coords[0] + _coords[3];
1067 #endif
1068 
1069 #if LIBMESH_DIM == 3
1070  return _coords[0] + _coords[4] + _coords[8];
1071 #endif
1072 }
template<typename T >
TypeTensor< T > libMesh::TypeTensor< T >::transpose ( ) const
inlineinherited

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

Definition at line 932 of file type_tensor.h.

933 {
934 #if LIBMESH_DIM == 1
935  return TypeTensor(_coords[0]);
936 #endif
937 
938 #if LIBMESH_DIM == 2
939  return TypeTensor(_coords[0],
940  _coords[2],
941  _coords[1],
942  _coords[3]);
943 #endif
944 
945 #if LIBMESH_DIM == 3
946  return TypeTensor(_coords[0],
947  _coords[3],
948  _coords[6],
949  _coords[1],
950  _coords[4],
951  _coords[7],
952  _coords[2],
953  _coords[5],
954  _coords[8]);
955 #endif
956 }
template<typename T >
void libMesh::TypeTensor< T >::write_unformatted ( std::ostream &  out,
const bool  newline = true 
) const
inherited

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

Definition at line 83 of file type_tensor.C.

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

Zero the tensor in any dimension.

Definition at line 1076 of file type_tensor.h.

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

1077 {
1078  for (unsigned int i=0; i<LIBMESH_DIM*LIBMESH_DIM; i++)
1079  _coords[i] = 0.;
1080 }

Member Data Documentation


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

Site Created By: libMesh Developers
Last modified: February 07 2014 16:58:00 UTC

Hosted By:
SourceForge.net Logo