www.mooseframework.org
Enumerations | Functions
RankTwoScalarTools Namespace Reference

Enumerations

enum  InvariantType {
  InvariantType::VonMisesStress, InvariantType::EffectiveStrain, InvariantType::Hydrostatic, InvariantType::L2norm,
  InvariantType::VolumetricStrain, InvariantType::FirstInvariant, InvariantType::SecondInvariant, InvariantType::ThirdInvariant,
  InvariantType::TriaxialityStress, InvariantType::MaxShear, InvariantType::StressIntensity, InvariantType::MaxPrincipal,
  InvariantType::MidPrincipal, InvariantType::MinPrincipal
}
 
enum  CylindricalComponent { CylindricalComponent::AxialStress, CylindricalComponent::HoopStress, CylindricalComponent::RadialStress }
 
enum  SphericalComponent { SphericalComponent::HoopStress, SphericalComponent::RadialStress }
 

Functions

MooseEnum scalarOptions ()
 This enum is left for legacy calls. More...
 
MooseEnum invariantOptions ()
 
MooseEnum cylindricalOptions ()
 
MooseEnum sphericalOptions ()
 
template<typename T >
component (const RankTwoTensorTempl< T > &r2tensor, unsigned int i, unsigned int j)
 
template<typename T >
component (const RankTwoTensorTempl< T > &r2tensor, unsigned int i, unsigned int j, Point &direction)
 
template<typename T >
vonMisesStress (const RankTwoTensorTempl< T > &stress)
 
template<typename T >
effectiveStrain (const RankTwoTensorTempl< T > &strain)
 
template<typename T >
hydrostatic (const RankTwoTensorTempl< T > &r2tensor)
 
template<typename T >
L2norm (const RankTwoTensorTempl< T > &r2tensor)
 
template<typename T >
volumetricStrain (const RankTwoTensorTempl< T > &strain)
 
template<typename T >
firstInvariant (const RankTwoTensorTempl< T > &r2tensor)
 
template<typename T >
secondInvariant (const RankTwoTensorTempl< T > &r2tensor)
 
template<typename T >
thirdInvariant (const RankTwoTensorTempl< T > &r2tensor)
 
template<typename T >
calcEigenValuesEigenVectors (const RankTwoTensorTempl< T > &r2tensor, unsigned int index, Point &eigenvec)
 
template<typename T >
maxPrincipal (const RankTwoTensorTempl< T > &r2tensor, Point &direction)
 
template<typename T >
midPrincipal (const RankTwoTensorTempl< T > &r2tensor, Point &direction)
 
template<typename T >
minPrincipal (const RankTwoTensorTempl< T > &r2tensor, Point &direction)
 
template<typename T >
axialStress (const RankTwoTensorTempl< T > &stress, const Point &point1, const Point &point2, Point &direction)
 
void normalPositionVector (const Point &point1, const Point &point2, const Point &curr_point, Point &normalPosition)
 
template<typename T >
hoopStress (const RankTwoTensorTempl< T > &stress, const Point &point1, const Point &point2, const Point &curr_point, Point &direction)
 
template<typename T >
radialStress (const RankTwoTensorTempl< T > &stress, const Point &point1, const Point &point2, const Point &curr_point, Point &direction)
 
template<typename T >
radialStress (const RankTwoTensorTempl< T > &stress, const Point &center, const Point &curr_point, Point &direction)
 
template<typename T >
hoopStress (const RankTwoTensorTempl< T > &stress, const Point &center, const Point &curr_point, Point &direction)
 
template<typename T >
directionValueTensor (const RankTwoTensorTempl< T > &r2tensor, const Point &direction)
 
template<typename T >
triaxialityStress (const RankTwoTensorTempl< T > &stress)
 
template<typename T >
maxShear (const RankTwoTensorTempl< T > &stress)
 
template<typename T >
stressIntensity (const RankTwoTensorTempl< T > &stress)
 
template<typename T >
getQuantity (const RankTwoTensorTempl< T > &tensor, const MooseEnum &scalar_type, const Point &point1, const Point &point2, const Point &curr_point, Point &direction)
 
template<typename T >
getCylindricalComponent (const RankTwoTensorTempl< T > &tensor, const CylindricalComponent &scalar_type, const Point &point1, const Point &point2, const Point &curr_point, Point &direction)
 
template<typename T >
getSphericalComponent (const RankTwoTensorTempl< T > &tensor, const SphericalComponent &scalar_type, const Point &center, const Point &curr_point, Point &direction)
 
template<typename T >
getPrincipalComponent (const RankTwoTensorTempl< T > &tensor, const InvariantType &scalar_type, Point &direction)
 
template<typename T >
getDirectionalComponent (const RankTwoTensorTempl< T > &tensor, const Point &direction)
 
template<typename T >
getInvariant (const RankTwoTensorTempl< T > &tensor, const InvariantType &scalar_type)
 
void setRotationMatrix (const RealVectorValue &outwardnormal, const RealVectorValue &axialVector, RankTwoTensor &rotationMatrix, const bool transpose)
 
template<bool is_ad>
void RankTwoTensorToVoigtVector (const GenericRankTwoTensor< is_ad > &tensor, GenericDenseVector< is_ad > &voigt_vector)
 
template<bool is_ad>
void VoigtVectorToRankTwoTensor (const GenericDenseVector< is_ad > &voigt_vector, GenericRankTwoTensor< is_ad > &tensor)
 

Enumeration Type Documentation

◆ CylindricalComponent

◆ InvariantType

◆ SphericalComponent

Function Documentation

◆ axialStress()

template<typename T >
T RankTwoScalarTools::axialStress ( const RankTwoTensorTempl< T > &  stress,
const Point &  point1,
const Point &  point2,
Point &  direction 
)

Definition at line 294 of file RankTwoScalarTools.h.

Referenced by RankTwoBasedFailureCriteriaNOSPD::computeFailureCriterionValue(), getCylindricalComponent(), and getQuantity().

298 {
299  Point axis = point2 - point1;
300  axis /= axis.norm();
301 
302  // Calculate the stress in the direction of the axis specifed by the user
303  T axial_stress = axis * (stress * axis);
304 
305  direction = axis;
306 
307  return axial_stress;
308 }
static const std::string axis
Definition: NS.h:27

◆ calcEigenValuesEigenVectors()

template<typename T >
T RankTwoScalarTools::calcEigenValuesEigenVectors ( const RankTwoTensorTempl< T > &  r2tensor,
unsigned int  index,
Point &  eigenvec 
)

Definition at line 231 of file RankTwoScalarTools.h.

Referenced by maxPrincipal(), midPrincipal(), minPrincipal(), and TEST().

234 {
235  std::vector<T> eigenval(LIBMESH_DIM);
236  RankTwoTensorTempl<T> eigvecs;
237  r2tensor.symmetricEigenvaluesEigenvectors(eigenval, eigvecs);
238 
239  T val = eigenval[index];
240  eigenvec = eigvecs.column(index);
241 
242  return val;
243 }
VectorValue< T > column(const unsigned int i) const
void symmetricEigenvaluesEigenvectors(std::vector< T > &eigvals, RankTwoTensorTempl< T > &eigvecs) const

◆ component() [1/2]

template<typename T >
T RankTwoScalarTools::component ( const RankTwoTensorTempl< T > &  r2tensor,
unsigned int  i,
unsigned int  j 
)

Definition at line 67 of file RankTwoScalarTools.h.

Referenced by MaterialTensorIntegralTempl< is_ad >::computeQpIntegral(), RankTwoCartesianComponentTempl< is_ad >::computeQpProperties(), RankTwoAuxTempl< is_ad >::computeValue(), NodalRankTwoPD::computeValue(), and LineMaterialRankTwoSampler::getScalarFromProperty().

68 {
69  return r2tensor(i, j);
70 }
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")

◆ component() [2/2]

template<typename T >
T RankTwoScalarTools::component ( const RankTwoTensorTempl< T > &  r2tensor,
unsigned int  i,
unsigned int  j,
Point &  direction 
)

Definition at line 78 of file RankTwoScalarTools.h.

79 {
80  direction.zero();
81  if (i == j)
82  direction(i) = 1.0;
83  else
84  {
85  direction(i) = std::sqrt(0.5);
86  direction(j) = std::sqrt(0.5);
87  }
88 
89  return r2tensor(i, j);
90 }
ADRealEigenVector< T, D, asd > sqrt(const ADRealEigenVector< T, D, asd > &)
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")

◆ cylindricalOptions()

MooseEnum RankTwoScalarTools::cylindricalOptions ( )

Definition at line 34 of file RankTwoScalarTools.C.

35 {
36  return MooseEnum("AxialStress HoopStress RadialStress");
37 }

◆ directionValueTensor()

template<typename T >
T RankTwoScalarTools::directionValueTensor ( const RankTwoTensorTempl< T > &  r2tensor,
const Point &  direction 
)

Definition at line 476 of file RankTwoScalarTools.h.

Referenced by getDirectionalComponent(), and getQuantity().

477 {
478  T tensor_value_in_direction = direction * (r2tensor * direction);
479 
480  return tensor_value_in_direction;
481 }

◆ effectiveStrain()

template<typename T >
T RankTwoScalarTools::effectiveStrain ( const RankTwoTensorTempl< T > &  strain)

Definition at line 113 of file RankTwoScalarTools.h.

Referenced by RankTwoBasedFailureCriteriaNOSPD::computeFailureCriterionValue(), and RankTwoInvariantTempl< is_ad >::computeQpProperties().

114 {
115  return std::sqrt(2.0 / 3.0 * strain.doubleContraction(strain));
116 }
ADRealEigenVector< T, D, asd > sqrt(const ADRealEigenVector< T, D, asd > &)
T doubleContraction(const RankTwoTensorTempl< T > &a) const

◆ firstInvariant()

template<typename T >
T RankTwoScalarTools::firstInvariant ( const RankTwoTensorTempl< T > &  r2tensor)

Definition at line 175 of file RankTwoScalarTools.h.

Referenced by getInvariant(), and getQuantity().

176 {
177  return r2tensor.trace();
178 }

◆ getCylindricalComponent()

template<typename T >
T RankTwoScalarTools::getCylindricalComponent ( const RankTwoTensorTempl< T > &  tensor,
const CylindricalComponent scalar_type,
const Point &  point1,
const Point &  point2,
const Point &  curr_point,
Point &  direction 
)

Definition at line 580 of file RankTwoScalarTools.h.

Referenced by RankTwoCylindricalComponentTempl< is_ad >::computeQpProperties().

586 {
587  switch (scalar_type)
588  {
589  case CylindricalComponent::AxialStress:
590  return axialStress(tensor, point1, point2, direction);
591  case CylindricalComponent::HoopStress:
592  return hoopStress(tensor, point1, point2, curr_point, direction);
593  case CylindricalComponent::RadialStress:
594  return radialStress(tensor, point1, point2, curr_point, direction);
595  default:
596  mooseError("RankTwoCylindricalComponent Error: scalar type invalid");
597  }
598 }
T axialStress(const RankTwoTensorTempl< T > &stress, const Point &point1, const Point &point2, Point &direction)
void mooseError(Args &&... args)
T hoopStress(const RankTwoTensorTempl< T > &stress, const Point &center, const Point &curr_point, Point &direction)
T radialStress(const RankTwoTensorTempl< T > &stress, const Point &center, const Point &curr_point, Point &direction)

◆ getDirectionalComponent()

template<typename T >
T RankTwoScalarTools::getDirectionalComponent ( const RankTwoTensorTempl< T > &  tensor,
const Point &  direction 
)

Definition at line 640 of file RankTwoScalarTools.h.

Referenced by RankTwoDirectionalComponentTempl< is_ad >::computeQpProperties().

641 {
642  return directionValueTensor(tensor, direction);
643 }
T directionValueTensor(const RankTwoTensorTempl< T > &r2tensor, const Point &direction)

◆ getInvariant()

template<typename T >
T RankTwoScalarTools::getInvariant ( const RankTwoTensorTempl< T > &  tensor,
const InvariantType scalar_type 
)

Definition at line 647 of file RankTwoScalarTools.h.

Referenced by RankTwoInvariantTempl< is_ad >::computeQpProperties().

648 {
649  switch (scalar_type)
650  {
651  case InvariantType::VonMisesStress:
652  return vonMisesStress(tensor);
653  case InvariantType::EffectiveStrain:
654  mooseError("To compute an effective inelastic strain use "
655  "RankTwoScalarTools::effectiveStrain()");
656  case InvariantType::Hydrostatic:
657  return hydrostatic(tensor);
659  return L2norm(tensor);
660  case InvariantType::VolumetricStrain:
661  return volumetricStrain(tensor);
662  case InvariantType::FirstInvariant:
663  return firstInvariant(tensor);
664  case InvariantType::SecondInvariant:
665  return secondInvariant(tensor);
666  case InvariantType::ThirdInvariant:
667  return thirdInvariant(tensor);
668  case InvariantType::TriaxialityStress:
669  return triaxialityStress(tensor);
670  case InvariantType::MaxShear:
671  return maxShear(tensor);
672  case InvariantType::StressIntensity:
673  return stressIntensity(tensor);
674  default:
675  mooseError("RankTwoCartesianComponent Error: Pass valid invariant - " +
676  invariantOptions().getRawNames());
677  }
678 }
MooseEnum invariantOptions()
T L2norm(const RankTwoTensorTempl< T > &r2tensor)
void mooseError(Args &&... args)
T thirdInvariant(const RankTwoTensorTempl< T > &r2tensor)
T hydrostatic(const RankTwoTensorTempl< T > &r2tensor)
T firstInvariant(const RankTwoTensorTempl< T > &r2tensor)
T vonMisesStress(const RankTwoTensorTempl< T > &stress)
T stressIntensity(const RankTwoTensorTempl< T > &stress)
T volumetricStrain(const RankTwoTensorTempl< T > &strain)
T maxShear(const RankTwoTensorTempl< T > &stress)
T secondInvariant(const RankTwoTensorTempl< T > &r2tensor)
T triaxialityStress(const RankTwoTensorTempl< T > &stress)

◆ getPrincipalComponent()

template<typename T >
T RankTwoScalarTools::getPrincipalComponent ( const RankTwoTensorTempl< T > &  tensor,
const InvariantType scalar_type,
Point &  direction 
)

Definition at line 621 of file RankTwoScalarTools.h.

Referenced by RankTwoInvariantTempl< is_ad >::computeQpProperties().

624 {
625  switch (scalar_type)
626  {
627  case InvariantType::MaxPrincipal:
628  return maxPrincipal(tensor, direction);
629  case InvariantType::MidPrincipal:
630  return midPrincipal(tensor, direction);
631  case InvariantType::MinPrincipal:
632  return minPrincipal(tensor, direction);
633  default:
634  mooseError("RankTwoInvariant Error: valid invariant");
635  }
636 }
T midPrincipal(const RankTwoTensorTempl< T > &r2tensor, Point &direction)
void mooseError(Args &&... args)
T maxPrincipal(const RankTwoTensorTempl< T > &r2tensor, Point &direction)
T minPrincipal(const RankTwoTensorTempl< T > &r2tensor, Point &direction)

◆ getQuantity()

template<typename T >
T RankTwoScalarTools::getQuantity ( const RankTwoTensorTempl< T > &  tensor,
const MooseEnum scalar_type,
const Point &  point1,
const Point &  point2,
const Point &  curr_point,
Point &  direction 
)

Definition at line 527 of file RankTwoScalarTools.h.

Referenced by RankTwoScalarAuxTempl< is_ad >::computeValue(), NodalRankTwoPD::computeValue(), XFEMRankTwoTensorMarkerUserObject::doesElementCrack(), MeshCut2DRankTwoTensorNucleation::doesElementCrack(), NodalRankTwoScalarPD::gatherWeightedValue(), and LineMaterialRankTwoScalarSampler::getScalarFromProperty().

533 {
534  switch (scalar_type)
535  {
536  case 0:
537  return vonMisesStress(tensor);
538  case 1:
539  mooseError("To compute an effective inelastic strain use "
540  "RankTwoScalarTools::effectiveStrain()");
541  case 2:
542  return hydrostatic(tensor);
543  case 3:
544  return L2norm(tensor);
545  case 4:
546  return maxPrincipal(tensor, direction);
547  case 5:
548  return midPrincipal(tensor, direction);
549  case 6:
550  return minPrincipal(tensor, direction);
551  case 7:
552  return volumetricStrain(tensor);
553  case 8:
554  return firstInvariant(tensor);
555  case 9:
556  return secondInvariant(tensor);
557  case 10:
558  return thirdInvariant(tensor);
559  case 11:
560  return axialStress(tensor, point1, point2, direction);
561  case 12:
562  return hoopStress(tensor, point1, point2, curr_point, direction);
563  case 13:
564  return radialStress(tensor, point1, point2, curr_point, direction);
565  case 14:
566  return triaxialityStress(tensor);
567  case 15:
568  return directionValueTensor(tensor, direction);
569  case 16:
570  return maxShear(tensor);
571  case 17:
572  return stressIntensity(tensor);
573  default:
574  mooseError("RankTwoScalarTools Error: invalid scalar type");
575  }
576 }
T axialStress(const RankTwoTensorTempl< T > &stress, const Point &point1, const Point &point2, Point &direction)
T midPrincipal(const RankTwoTensorTempl< T > &r2tensor, Point &direction)
T L2norm(const RankTwoTensorTempl< T > &r2tensor)
void mooseError(Args &&... args)
T hoopStress(const RankTwoTensorTempl< T > &stress, const Point &center, const Point &curr_point, Point &direction)
T maxPrincipal(const RankTwoTensorTempl< T > &r2tensor, Point &direction)
T thirdInvariant(const RankTwoTensorTempl< T > &r2tensor)
T hydrostatic(const RankTwoTensorTempl< T > &r2tensor)
T firstInvariant(const RankTwoTensorTempl< T > &r2tensor)
T vonMisesStress(const RankTwoTensorTempl< T > &stress)
T stressIntensity(const RankTwoTensorTempl< T > &stress)
T volumetricStrain(const RankTwoTensorTempl< T > &strain)
T radialStress(const RankTwoTensorTempl< T > &stress, const Point &center, const Point &curr_point, Point &direction)
T minPrincipal(const RankTwoTensorTempl< T > &r2tensor, Point &direction)
T maxShear(const RankTwoTensorTempl< T > &stress)
T secondInvariant(const RankTwoTensorTempl< T > &r2tensor)
T directionValueTensor(const RankTwoTensorTempl< T > &r2tensor, const Point &direction)
T triaxialityStress(const RankTwoTensorTempl< T > &stress)

◆ getSphericalComponent()

template<typename T >
T RankTwoScalarTools::getSphericalComponent ( const RankTwoTensorTempl< T > &  tensor,
const SphericalComponent scalar_type,
const Point &  center,
const Point &  curr_point,
Point &  direction 
)

Definition at line 602 of file RankTwoScalarTools.h.

Referenced by RankTwoSphericalComponentTempl< is_ad >::computeQpProperties().

607 {
608  switch (scalar_type)
609  {
610  case SphericalComponent::HoopStress:
611  return hoopStress(tensor, center, curr_point, direction);
612  case SphericalComponent::RadialStress:
613  return radialStress(tensor, center, curr_point, direction);
614  default:
615  mooseError("RankTwoSphericalComponent Error: scalar type invalid");
616  }
617 }
void mooseError(Args &&... args)
T hoopStress(const RankTwoTensorTempl< T > &stress, const Point &center, const Point &curr_point, Point &direction)
T radialStress(const RankTwoTensorTempl< T > &stress, const Point &center, const Point &curr_point, Point &direction)
static const std::string center
Definition: NS.h:28

◆ hoopStress() [1/2]

template<typename T >
T RankTwoScalarTools::hoopStress ( const RankTwoTensorTempl< T > &  stress,
const Point &  point1,
const Point &  point2,
const Point &  curr_point,
Point &  direction 
)

Definition at line 336 of file RankTwoScalarTools.h.

Referenced by getCylindricalComponent(), getQuantity(), and getSphericalComponent().

341 {
342  // Calculate the cross of the normal to the axis of rotation from the current
343  // location and the axis of rotation
344  Point xp;
345  normalPositionVector(point1, point2, curr_point, xp);
346  Point axis_rotation = point2 - point1;
347  Point yp = axis_rotation / axis_rotation.norm();
348  Point zp = xp.cross(yp);
349 
350  // Calculate the scalar value of the hoop stress
351  T hoop_stress = zp * (stress * zp);
352 
353  direction = zp;
354 
355  return hoop_stress;
356 }
void normalPositionVector(const Point &point1, const Point &point2, const Point &curr_point, Point &normalPosition)

◆ hoopStress() [2/2]

template<typename T >
T RankTwoScalarTools::hoopStress ( const RankTwoTensorTempl< T > &  stress,
const Point &  center,
const Point &  curr_point,
Point &  direction 
)

Definition at line 422 of file RankTwoScalarTools.h.

426 {
427  Point radial = curr_point - center;
428  Real r = radial.norm();
429  radial /= r;
430 
431  // Given normal vector N=(n1,n2,n3) and current point C(c1,c2,c3), the tangential plane is then
432  // defined as n1(x-c1) + n2(y-c2) + n3(z-c3)=0. Let us assume n1!=0, the arbitrary point P on this
433  // plane can be taken as P(x,c2+r,c3+r) where r is the radius. The x can be solved as x =
434  // -r(n2+n3)/n1 + c1. The tangential vector PC is given as P-C.
435 
436  Point tangential;
437 
438  if (!MooseUtils::absoluteFuzzyEqual(radial(0), 0.0))
439  {
440  Real x = curr_point(0) - (radial(1) + radial(2)) * r / radial(0);
441  tangential = Point(x, curr_point(1) + r, curr_point(2) + r) - curr_point;
442  }
443  else if (!MooseUtils::absoluteFuzzyEqual(radial(1), 0.0))
444  {
445  Real y = curr_point(1) - (radial(0) + radial(2)) * r / radial(1);
446  tangential = Point(curr_point(0) + r, y, curr_point(2) + r) - curr_point;
447  }
448  else if (!MooseUtils::absoluteFuzzyEqual(radial(1), 0.0))
449  {
450  Real z = curr_point(2) - (radial(0) + radial(1)) * r / radial(2);
451  tangential = Point(curr_point(0) + r, curr_point(1) + r, z) - curr_point;
452  }
453  else
454  mooseError("In Hoop stress calculation for spherical geometry, the current (quadracture) point "
455  "is likely to be at the center. ");
456 
457  tangential /= tangential.norm();
458 
459  // Compute the scalar stress component in the hoop direction
460  T hoop_stress = 0.0;
461  for (unsigned int i = 0; i < 3; ++i)
462  for (unsigned int j = 0; j < 3; ++j)
463  hoop_stress += tangential(j) * stress(j, i) * tangential(i);
464 
465  direction = tangential;
466 
467  return hoop_stress;
468 }
bool absoluteFuzzyEqual(const T &var1, const T2 &var2, const T3 &tol=libMesh::TOLERANCE *libMesh::TOLERANCE)
void mooseError(Args &&... args)
const std::vector< double > y
const std::vector< double > x
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
static const std::string center
Definition: NS.h:28

◆ hydrostatic()

template<typename T >
T RankTwoScalarTools::hydrostatic ( const RankTwoTensorTempl< T > &  r2tensor)

Definition at line 124 of file RankTwoScalarTools.h.

Referenced by getInvariant(), getQuantity(), and triaxialityStress().

125 {
126  return r2tensor.trace() / 3.0;
127 }

◆ invariantOptions()

MooseEnum RankTwoScalarTools::invariantOptions ( )

Definition at line 26 of file RankTwoScalarTools.C.

Referenced by getInvariant().

27 {
28  return MooseEnum("VonMisesStress EffectiveStrain Hydrostatic L2norm VolumetricStrain "
29  "FirstInvariant SecondInvariant "
30  "ThirdInvariant TriaxialityStress MaxShear StressIntensity EffectiveStrain");
31 }

◆ L2norm()

template<typename T >
T RankTwoScalarTools::L2norm ( const RankTwoTensorTempl< T > &  r2tensor)

◆ maxPrincipal()

template<typename T >
T RankTwoScalarTools::maxPrincipal ( const RankTwoTensorTempl< T > &  r2tensor,
Point &  direction 
)

Definition at line 253 of file RankTwoScalarTools.h.

Referenced by RankTwoBasedFailureCriteriaNOSPD::computeFailureCriterionValue(), getPrincipalComponent(), getQuantity(), and maxShear().

254 {
255  return calcEigenValuesEigenVectors(r2tensor, (LIBMESH_DIM - 1), direction);
256 }
T calcEigenValuesEigenVectors(const RankTwoTensorTempl< T > &r2tensor, unsigned int index, Point &eigenvec)

◆ maxShear()

template<typename T >
T RankTwoScalarTools::maxShear ( const RankTwoTensorTempl< T > &  stress)

Definition at line 499 of file RankTwoScalarTools.h.

Referenced by RankTwoBasedFailureCriteriaNOSPD::computeFailureCriterionValue(), getInvariant(), getQuantity(), and stressIntensity().

500 {
501  Point dummy;
502  return (maxPrincipal(stress, dummy) - minPrincipal(stress, dummy)) / 2.;
503 }
T maxPrincipal(const RankTwoTensorTempl< T > &r2tensor, Point &direction)
T minPrincipal(const RankTwoTensorTempl< T > &r2tensor, Point &direction)

◆ midPrincipal()

template<typename T >
T RankTwoScalarTools::midPrincipal ( const RankTwoTensorTempl< T > &  r2tensor,
Point &  direction 
)

Definition at line 266 of file RankTwoScalarTools.h.

Referenced by getPrincipalComponent(), and getQuantity().

267 {
268  return calcEigenValuesEigenVectors(r2tensor, 1, direction);
269 }
T calcEigenValuesEigenVectors(const RankTwoTensorTempl< T > &r2tensor, unsigned int index, Point &eigenvec)

◆ minPrincipal()

template<typename T >
T RankTwoScalarTools::minPrincipal ( const RankTwoTensorTempl< T > &  r2tensor,
Point &  direction 
)

Definition at line 279 of file RankTwoScalarTools.h.

Referenced by getPrincipalComponent(), getQuantity(), and maxShear().

280 {
281  return calcEigenValuesEigenVectors(r2tensor, 0, direction);
282 }
T calcEigenValuesEigenVectors(const RankTwoTensorTempl< T > &r2tensor, unsigned int index, Point &eigenvec)

◆ normalPositionVector()

void RankTwoScalarTools::normalPositionVector ( const Point &  point1,
const Point &  point2,
const Point &  curr_point,
Point &  normalPosition 
)

Definition at line 46 of file RankTwoScalarTools.C.

Referenced by hoopStress(), and radialStress().

50 {
51  // Find the nearest point on the axis of rotation (defined by point2 - point1)
52  // to the current position, e.g. the normal to the axis of rotation at the
53  // current position
54  Point axis_rotation = point2 - point1;
55  Point positionWRTpoint1 = point1 - curr_point;
56  Real projection = (axis_rotation * positionWRTpoint1) / axis_rotation.norm_sq();
57  Point normal = point1 - projection * axis_rotation;
58 
59  // Calculate the direction normal to the plane formed by the axis of rotation
60  // and the normal to the axis of rotation from the current position.
61  normalPosition = curr_point - normal;
62  normalPosition /= normalPosition.norm();
63 }
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real

◆ radialStress() [1/2]

template<typename T >
T RankTwoScalarTools::radialStress ( const RankTwoTensorTempl< T > &  stress,
const Point &  point1,
const Point &  point2,
const Point &  curr_point,
Point &  direction 
)

Definition at line 369 of file RankTwoScalarTools.h.

Referenced by getCylindricalComponent(), getQuantity(), and getSphericalComponent().

374 {
375  Point radial_norm;
376  normalPositionVector(point1, point2, curr_point, radial_norm);
377 
378  // Compute the scalar stress component in the direction of the normal vector from the
379  // user-defined axis of rotation.
380  T radial_stress = radial_norm * (stress * radial_norm);
381 
382  direction = radial_norm;
383 
384  return radial_stress;
385 }
void normalPositionVector(const Point &point1, const Point &point2, const Point &curr_point, Point &normalPosition)

◆ radialStress() [2/2]

template<typename T >
T RankTwoScalarTools::radialStress ( const RankTwoTensorTempl< T > &  stress,
const Point &  center,
const Point &  curr_point,
Point &  direction 
)

Definition at line 397 of file RankTwoScalarTools.h.

401 {
402  Point radial = curr_point - center;
403  radial /= radial.norm();
404 
405  // Compute the scalar stress component in the radial direction
406  T radial_stress = radial * (stress * radial);
407 
408  direction = radial;
409 
410  return radial_stress;
411 }
static const std::string center
Definition: NS.h:28

◆ RankTwoTensorToVoigtVector()

template<bool is_ad>
void RankTwoScalarTools::RankTwoTensorToVoigtVector ( const GenericRankTwoTensor< is_ad > &  tensor,
GenericDenseVector< is_ad > &  voigt_vector 
)

Definition at line 687 of file RankTwoScalarTools.h.

689 {
690  voigt_vector(0) = tensor(0, 0);
691  voigt_vector(1) = tensor(1, 1);
692  voigt_vector(2) = tensor(2, 2);
693  voigt_vector(3) = tensor(0, 1);
694  voigt_vector(4) = tensor(1, 2);
695  voigt_vector(5) = tensor(0, 2);
696 }

◆ scalarOptions()

MooseEnum RankTwoScalarTools::scalarOptions ( )

This enum is left for legacy calls.

Definition at line 17 of file RankTwoScalarTools.C.

Referenced by MeshCut2DRankTwoTensorNucleation::validParams(), XFEMRankTwoTensorMarkerUserObject::validParams(), NodalRankTwoScalarPD::validParams(), RankTwoScalarAuxTempl< is_ad >::validParams(), LineMaterialRankTwoScalarSampler::validParams(), and NodalRankTwoPD::validParams().

18 {
19  return MooseEnum("VonMisesStress EffectiveStrain Hydrostatic L2norm MaxPrincipal "
20  "MidPrincipal MinPrincipal VolumetricStrain FirstInvariant SecondInvariant "
21  "ThirdInvariant AxialStress HoopStress RadialStress TriaxialityStress "
22  "Direction MaxShear StressIntensity");
23 }

◆ secondInvariant()

template<typename T >
T RankTwoScalarTools::secondInvariant ( const RankTwoTensorTempl< T > &  r2tensor)

Definition at line 189 of file RankTwoScalarTools.h.

Referenced by getInvariant(), and getQuantity().

190 {
191  T val = 0.0;
192  for (unsigned int i = 0; i < 2; ++i)
193  for (unsigned int j = i + 1; j < 3; ++j)
194  {
195  val += r2tensor(i, i) * r2tensor(j, j);
196  val -= (r2tensor(i, j) * r2tensor(i, j) + r2tensor(j, i) * r2tensor(j, i)) * 0.5;
197  }
198 
199  return val;
200 }
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")

◆ setRotationMatrix()

void RankTwoScalarTools::setRotationMatrix ( const RealVectorValue outwardnormal,
const RealVectorValue axialVector,
RankTwoTensor rotationMatrix,
const bool  transpose 
)

Definition at line 66 of file RankTwoScalarTools.C.

Referenced by HillElastoPlasticityStressUpdateTempl< is_ad >::initQpStatefulProperties().

70 {
71  RealVectorValue radial_vector(outwardnormal);
72  RealVectorValue azimuthal_vector(0, 0, 0);
73  RealVectorValue polar_vector(axialVector);
74 
75  Real rv_norm = radial_vector.norm();
76  if (rv_norm)
77  radial_vector /= rv_norm;
78  else
79  mooseError("The outward normal vector cannot be a zero vector");
80 
81  azimuthal_vector = polar_vector.cross(radial_vector);
82  azimuthal_vector /= azimuthal_vector.norm();
83 
84  if (!transpose)
85  {
86  // Considering that Moose uses convention [T'] = [Q][T][Q_transpose], the
87  // basis vectors of old coordinate system should be arranged column wise
88  // to construct the rotation matrix Q.
89 
90  rotationMatrix.fillColumn(0, radial_vector);
91  rotationMatrix.fillColumn(1, azimuthal_vector);
92  rotationMatrix.fillColumn(2, polar_vector);
93  }
94  else
95  {
96  rotationMatrix.fillRow(0, radial_vector);
97  rotationMatrix.fillRow(1, azimuthal_vector);
98  rotationMatrix.fillRow(2, polar_vector);
99  }
100 }
void fillRow(unsigned int r, const libMesh::TypeVector< Real > &v)
void mooseError(Args &&... args)
void fillColumn(unsigned int c, const libMesh::TypeVector< Real > &v)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real

◆ sphericalOptions()

MooseEnum RankTwoScalarTools::sphericalOptions ( )

Definition at line 40 of file RankTwoScalarTools.C.

41 {
42  return MooseEnum("HoopStress RadialStress");
43 }

◆ stressIntensity()

template<typename T >
T RankTwoScalarTools::stressIntensity ( const RankTwoTensorTempl< T > &  stress)

Definition at line 509 of file RankTwoScalarTools.h.

Referenced by getInvariant(), and getQuantity().

510 {
511  return 2. * maxShear(stress);
512 }
T maxShear(const RankTwoTensorTempl< T > &stress)

◆ thirdInvariant()

template<typename T >
T RankTwoScalarTools::thirdInvariant ( const RankTwoTensorTempl< T > &  r2tensor)

Definition at line 208 of file RankTwoScalarTools.h.

Referenced by getInvariant(), and getQuantity().

209 {
210  T val = 0.0;
211  val = r2tensor(0, 0) * r2tensor(1, 1) * r2tensor(2, 2) -
212  r2tensor(0, 0) * r2tensor(1, 2) * r2tensor(2, 1) +
213  r2tensor(0, 1) * r2tensor(1, 2) * r2tensor(2, 0) -
214  r2tensor(0, 1) * r2tensor(1, 0) * r2tensor(2, 2) +
215  r2tensor(0, 2) * r2tensor(1, 0) * r2tensor(2, 1) -
216  r2tensor(0, 2) * r2tensor(1, 1) * r2tensor(2, 0);
217 
218  return val;
219 }

◆ triaxialityStress()

template<typename T >
T RankTwoScalarTools::triaxialityStress ( const RankTwoTensorTempl< T > &  stress)

Definition at line 488 of file RankTwoScalarTools.h.

Referenced by getInvariant(), and getQuantity().

489 {
490  return hydrostatic(stress) / vonMisesStress(stress);
491 }
T hydrostatic(const RankTwoTensorTempl< T > &r2tensor)
T vonMisesStress(const RankTwoTensorTempl< T > &stress)

◆ VoigtVectorToRankTwoTensor()

template<bool is_ad>
void RankTwoScalarTools::VoigtVectorToRankTwoTensor ( const GenericDenseVector< is_ad > &  voigt_vector,
GenericRankTwoTensor< is_ad > &  tensor 
)

Definition at line 700 of file RankTwoScalarTools.h.

702 {
703  tensor(0, 0) = voigt_vector(0);
704  tensor(1, 1) = voigt_vector(1);
705  tensor(2, 2) = voigt_vector(2);
706  tensor(0, 1) = tensor(1, 0) = voigt_vector(3);
707  tensor(1, 2) = tensor(2, 1) = voigt_vector(4);
708  tensor(0, 2) = tensor(2, 0) = voigt_vector(5);
709 }

◆ volumetricStrain()

template<typename T >
T RankTwoScalarTools::volumetricStrain ( const RankTwoTensorTempl< T > &  strain)

Definition at line 164 of file RankTwoScalarTools.h.

Referenced by getInvariant(), and getQuantity().

165 {
166  return std::exp(strain(0, 0)) * std::exp(strain(1, 1)) * std::exp(strain(2, 2)) - 1.0;
167 }

◆ vonMisesStress()

template<typename T >
T RankTwoScalarTools::vonMisesStress ( const RankTwoTensorTempl< T > &  stress)

Definition at line 100 of file RankTwoScalarTools.h.

Referenced by RankTwoBasedFailureCriteriaNOSPD::computeFailureCriterionValue(), getInvariant(), getQuantity(), and triaxialityStress().

101 {
102  RankTwoTensorTempl<T> dev_stress = stress.deviatoric();
103 
104  return std::sqrt(1.5 * dev_stress.doubleContraction(dev_stress));
105 }
ADRealEigenVector< T, D, asd > sqrt(const ADRealEigenVector< T, D, asd > &)
RankTwoTensorTempl< T > deviatoric() const
T doubleContraction(const RankTwoTensorTempl< T > &a) const