22 #include "libmesh/libmesh.h" 23 #include "libmesh/tensor_value.h" 24 #include "libmesh/vector_value.h" 25 #include "libmesh/int_range.h" 27 #include "metaphysicl/raw_type.h" 92 static constexpr
unsigned int N2 =
N *
N;
138 void print(std::ostream & stm = Moose::out)
const;
141 void printReal(std::ostream & stm = Moose::out)
const;
144 void printDualReal(
unsigned int nDual, std::ostream & stm = Moose::out)
const;
224 const T & S11,
const T & S22,
const T & S33,
const T & S23,
const T & S13,
const T & S12);
272 template <
typename T2>
287 template <
typename T2>
476 VectorValue<T>
column(
const unsigned int i)
const;
543 template <
typename Scalar>
544 typename boostcopy::enable_if_c<ScalarTraits<Scalar>::value,
RankTwoTensorTempl &>::type
547 libmesh_assert_equal_to(p, Scalar(0));
836 void addIa(
const T & a);
942 template <
typename T2>
964 template <
typename T2>
998 template <typename T2, typename std::enable_if<ScalarTraits<T2>::value,
int>::type = 0>
1015 template <typename T2, typename std::enable_if<ScalarTraits<T2>::value,
int>::type = 0>
1036 template <
typename T2>
1058 template <
typename T2>
1114 template <
int n,
int o,
int p,
int q>
1144 template <
int n,
int o,
int p,
int q,
int r,
int s>
1152 usingTensorIndices(i_, j_, k_, l_);
1153 return times<i_, j_, k_, l_>(b);
1246 T
sin3Lode(
const T & r0,
const T & r0_value)
const;
1358 void syev(
const char * calculation_type, std::vector<T> & eigvals, std::vector<T> & a)
const;
1378 template <
typename T>
1395 template <
typename T>
1396 template <
typename T2>
1403 template <
typename T>
1404 template <
typename T2>
1411 template <
typename T>
1412 template <typename T2, typename std::enable_if<ScalarTraits<T2>::value,
int>::type>
1419 template <
typename T>
1420 template <
typename T2>
1427 template <
typename T>
1428 template <
typename T2>
1435 template <
typename T>
1436 template <typename T2, typename std::enable_if<ScalarTraits<T2>::value,
int>::type>
1443 template <
typename T>
1451 this->symmetricEigenvaluesEigenvectors(eigval, eigvec);
1454 std::array<T, N> epos;
1458 epos[i] = (
std::abs(eigval[i]) + eigval[i]) / 2.0;
1459 d[i] = 0 < eigval[i] ? 1.0 : 0.0;
1469 proj_pos += d[a] * Ma.outerProduct(Ma);
1472 usingTensorIndices(i_, j_, k_, l_);
1479 Gab = Ma.template times<i_, k_, j_, l_>(Mb) + Ma.template times<i_, l_, j_, k_>(Mb);
1480 Gba = Mb.template times<i_, k_, j_, l_>(Ma) + Mb.template times<i_, l_, j_, k_>(Ma);
1484 theta_ab = 0.5 * (epos[a] - epos[b]) / (eigval[a] - eigval[b]);
1486 theta_ab = 0.25 * (d[a] + d[b]);
1488 proj_pos += theta_ab * (Gab + Gba);
1493 mooseError(
"positiveProjectionEigenDecomposition is only available for ordered tensor " 1497 template <
typename T>
1503 T bar = secondInvariant();
1513 mooseError(
"sin3Lode is only available for ordered tensor component types");
1516 template <
typename T>
1522 T bar = secondInvariant();
1527 (dthirdInvariant() /
std::pow(bar, 1.5) -
1528 1.5 * dsecondInvariant() * thirdInvariant() /
std::pow(bar, 2.5));
1531 mooseError(
"dsin3Lode is only available for ordered tensor component types");
1534 template <
typename T>
1540 T bar = secondInvariant();
1544 T J3 = thirdInvariant();
1548 1.5 * d2secondInvariant() * J3 /
std::pow(bar, 2.5);
1550 for (
unsigned i = 0; i <
N; ++i)
1551 for (
unsigned j = 0; j <
N; ++j)
1552 for (
unsigned k = 0; k <
N; ++k)
1553 for (
unsigned l = 0; l <
N; ++l)
1554 deriv(i, j, k, l) += (-1.5 * dII(i, j) * dIII(k, l) - 1.5 * dIII(i, j) * dII(k, l)) /
1556 1.5 * 2.5 * dII(i, j) * dII(k, l) * J3 /
std::pow(bar, 3.5);
1562 mooseError(
"d2sin3Lode is only available for ordered tensor component types");
1565 template <
typename T>
1566 template <
int n,
int o,
int p,
int q>
1572 for (x[0] = 0; x[0] <
N; ++x[0])
1573 for (x[1] = 0; x[1] <
N; ++x[1])
1574 for (x[2] = 0; x[2] <
N; ++x[2])
1575 for (x[3] = 0; x[3] <
N; ++x[3])
1576 result(x[0], x[1], x[2], x[3]) = (*this)(x[n], x[o]) * b(x[p], x[q]);
1581 template <
typename T>
1582 template <
int n,
int o,
int p,
int q,
int r,
int s>
1588 for (x[0] = 0; x[0] <
N; ++x[0])
1589 for (x[1] = 0; x[1] <
N; ++x[1])
1590 for (x[2] = 0; x[2] <
N; ++x[2])
1591 for (x[3] = 0; x[3] <
N; ++x[3])
1592 for (x[4] = 0; x[4] <
N; ++x[4])
1593 result(x[0], x[1], x[2], x[3]) += (*this)(x[n], x[o]) * b(x[p], x[q], x[r], x[s]);
RankFourTensorTempl< T > outerProduct(const RankTwoTensorTempl< T > &b) const
Return the outer product .
RankFourTensorTempl is designed to handle any N-dimensional fourth order tensor, C.
void fillRow(unsigned int r, const libMesh::TypeVector< T > &v)
Assign values to a specific row of the second order tensor.
RankTwoTensorTempl< T > dsecondInvariant() const
Return the derivative of the main second invariant w.r.t.
RankTwoTensorTempl< T > inverse() const
Return the inverse of this second order tensor.
RankFourTensorTempl< T > d2secondInvariant() const
Return the second derivative of the main second invariant w.r.t.
RankTwoTensorTempl(const RankTwoTensorTempl< T2 > &a)
The conversion operator from RankTwoTensorTempl<T2> to RankTwoTensorTempl<T> where T2 is convertible ...
RankTwoTensorTempl()
Empty constructor; fills to zero.
bool operator==(const RankTwoTensorTempl< T > &a) const
Defines logical equality with another RankTwoTensorTempl<T>
void mooseSetToZero(T &v)
Helper function templates to set a variable to zero.
T generalSecondInvariant() const
Return the principal second invariant of this second order tensor.
RankFourTensorTempl< T > positiveProjectionEigenDecomposition(std::vector< T > &, RankTwoTensorTempl< T > &) const
Return the positive projection tensor.
static RankTwoTensorTempl initializeSymmetric(const libMesh::TypeVector< T > &v0, const libMesh::TypeVector< T > &v1, const libMesh::TypeVector< T > &v2)
Named constructor for initializing symmetrically.
static void initRandom(unsigned int)
Initialize the random seed based on an unsigned integer.
void symmetricEigenvalues(std::vector< T > &eigvals) const
computes eigenvalues, assuming tens is symmetric, and places them in ascending order in eigvals ...
T sin3Lode(const T &r0, const T &r0_value) const
Sin(3*Lode_angle)
bool absoluteFuzzyEqual(const T &var1, const T2 &var2, const T3 &tol=libMesh::TOLERANCE *libMesh::TOLERANCE)
Function to check whether two variables are equal within an absolute tolerance.
This class defines a Tensor that can change its shape.
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
void dsymmetricEigenvalues(std::vector< T > &eigvals, std::vector< RankTwoTensorTempl< T >> &deigvals) const
computes eigenvalues, and their symmetric derivatives wrt vals, assuming tens is symmetric ...
static RankTwoTensorTempl< T > genRandomTensor(T stddev, T mean)
Generate a random second order tensor with all 9 components treated as independent random variables f...
void print(std::ostream &stm=Moose::out) const
Print the rank two tensor.
RankTwoTensorTempl< T > & operator-=(const RankTwoTensorTempl< T > &a)
Subtract another second order tensor from this one.
static RankTwoTensorTempl< T > plusTranspose(const RankTwoTensorTempl< T > &)
Initialize a second order tensor with expression .
RankThreeTensorTempl< T > mixedProductJkI(const VectorValue< T > &b) const
Return the tensor product of this second order tensor with a vector .
RankThreeTensorTempl< T > contraction(const RankThreeTensorTempl< T > &b) const
Return the single contraction of this second order tensor with a third order tensor ...
auto operator*(const Scalar &scalar) const -> typename boostcopy::enable_if_c< ScalarTraits< Scalar >::value, TypeTensor< decltype(T() *scalar)>>::type
VectorValue< T > column(const unsigned int i) const
Get the i-th column of the second order tensor.
void surfaceFillFromInputVector(const std::vector< T > &input)
sets _coords[0][0], _coords[0][1], _coords[1][0], _coords[1][1] to input, and the remainder to zero ...
void mooseSetToZero< ADRankTwoTensor >(ADRankTwoTensor &v)
Helper function template specialization to set an object to zero.
static constexpr std::size_t dim
This is the dimension of all vector and tensor datastructures used in MOOSE.
RankTwoTensorTempl< T > & operator+=(const RankTwoTensorTempl< T > &a)
Add another second order tensor to this one.
T secondInvariant() const
Return the main second invariant of this second order tensor.
RankTwoTensorTempl< T > dtrace() const
Return the derivative of the trace w.r.t.
ADRealEigenVector< T, D, asd > sqrt(const ADRealEigenVector< T, D, asd > &)
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
static RankTwoTensorTempl Identity()
Initialize a second order identity tensor.
RankTwoTensorTempl< T > dthirdInvariant() const
Denote the _coords[i][j] by A_ij, then this returns d(thirdInvariant()/dA_ij.
boostcopy::enable_if_c< ScalarTraits< Scalar >::value, TypeTensor< typename CompareTypes< T, Scalar >::supertype > >::type operator/(const Scalar &) const
static constexpr unsigned int N
Tensor dimension, i.e.
auto max(const L &left, const R &right)
void rotate(const RankTwoTensorTempl< T > &R)
Rotate the tensor in-place given a rotation tensor .
RankThreeTensor is designed to handle any N-dimensional third order tensor, r.
static constexpr unsigned int N2
The square of the tensor dimension.
static constexpr Real identityCoords[N2]
friend void dataLoad(std::istream &, RankTwoTensorTempl< T2 > &, void *)
T trace() const
A wrapper for tr()
void fillColumn(unsigned int c, const libMesh::TypeVector< T > &v)
Assign values to a specific column of the second order tensor.
friend void dataStore(std::ostream &, RankTwoTensorTempl< T2 > &, void *)
ADRealEigenVector< T, D, asd > abs(const ADRealEigenVector< T, D, asd > &)
RankTwoTensorTempl< T > dsin3Lode(const T &r0) const
d(sin3Lode)/dA_ij If secondInvariant() <= r0 then return zero This is to gaurd against precision-loss...
TypeTensor< typename CompareTypes< T, T2 >::supertype > operator+(const TypeTensor< T2 > &) const
static RankTwoTensorTempl< T > selfOuterProduct(const libMesh::TypeVector< T > &)
Initialize a second order tensor as the outer product of a vector with itself, i.e.
static RankTwoTensorTempl initializeFromRows(const libMesh::TypeVector< T > &row0, const libMesh::TypeVector< T > &row1, const libMesh::TypeVector< T > &row2)
Named constructor for initializing from row vectors.
RankTwoTensorTempl< T > deviatoric() const
Return the deviatoric part of this tensor .
static MooseEnum fillMethodEnum()
Get the available FillMethod options.
T L2norm() const
Sqrt(_coords[i][j]*_coords[i][j])
void fillFromInputVector(const std::vector< T > &input, FillMethod fill_method=autodetect)
The smart mutator that determines how to fill the second order tensor based on the size of the input ...
void setToIdentity()
Set the tensor to identity.
Custom type trait that has a value of true for types that cam be use interchangably with Real...
TypeTensor< T > operator-() const
RankTwoTensorTempl< typename CompareTypes< T, T2 >::supertype > operator*(const T2 &a) const
Return this tensor multiplied by a scalar (component-wise)
RankTwoTensorTempl< T > & operator/=(const T &a)
Divide this tensor by a scalar (component-wise)
Real deriv(unsigned n, unsigned alpha, unsigned beta, Real x)
RankFourTensorTempl< T > d2thirdInvariant() const
Denote the _coords[i][j] by A_ij, then this returns d^2(thirdInvariant)/dA_ij/dA_kl.
void addIa(const T &a)
Add identity times a to _coords.
RankTwoTensorTempl< T > rotateXyPlane(T a)
Rotate the tensor about the z-axis.
SymmetricRankTwoTensorTempl is designed to handle the Stress or Strain Tensor for an anisotropic mate...
void d2symmetricEigenvalues(std::vector< RankFourTensorTempl< T >> &deriv) const
Computes second derivatives of Eigenvalues of a rank two tensor.
static RankTwoTensorTempl< T > genRandomSymmTensor(T stddev, T mean)
Generate a random symmetric second order tensor with the 6 upper-triangular components treated as ind...
RankTwoTensorTempl< T > rotated(const RankTwoTensorTempl< T > &R) const
Return the rotated tensor given a rotation tensor .
T doubleContraction(const RankTwoTensorTempl< T > &a) const
Return the double contraction with another second order tensor .
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
RankTwoTensorTempl< T > square() const
Return .
RankFourTensorTempl< T > d2sin3Lode(const T &r0) const
d^2(sin3Lode)/dA_ij/dA_kl If secondInvariant() <= r0 then return zero This is to gaurd against precis...
RankTwoTensorTempl< T > & operator*=(const T &a)
Multiply this tensor by a scalar (component-wise)
T thirdInvariant() const
Denote the _coords[i][j] by A_ij, then S_ij = A_ij - de_ij*tr(A)/3 Then this returns det(S + S...
boostcopy::enable_if_c< ScalarTraits< Scalar >::value, RankTwoTensorTempl & >::type operator=(const Scalar &libmesh_dbg_var(p))
Assignment-from-scalar operator.
void printReal(std::ostream &stm=Moose::out) const
Print the Real part of the RankTwoTensorTempl<ADReal>
static RankTwoTensorTempl< T > timesTranspose(const RankTwoTensorTempl< T > &)
Initialize a second order tensor with expression .
OutputTools< Real >::VariableValue VariableValue
RankTwoTensorTempl< typename CompareTypes< T, T2 >::supertype > operator+(const libMesh::TypeTensor< T2 > &a) const
Return the sum of two second order tensors.
static RankTwoTensorTempl< T > transposeTimes(const RankTwoTensorTempl< T > &)
Initialize a second order tensor with expression .
RankTwoTensorTempl< T > transpose() const
Return the tensor transposed.
RankTwoTensorTempl< T > & operator=(const RankTwoTensorTempl< T > &a)
Assignment operator.
void vectorOuterProduct(const libMesh::TypeVector< T > &, const libMesh::TypeVector< T > &)
Set the values of the second order tensor to be the outer product of two vectors, i...
RankTwoTensorTempl(const std::vector< T > &input)
Constructor that proxies the fillFromInputVector() method.
void fillRealTensor(libMesh::TensorValue< T > &)
Fill a libMesh::TensorValue<T> from this second order tensor.
RankTwoTensorTempl(const libMesh::TypeTensor< T > &a)
The conversion operator from a libMesh::TypeTensor
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
bool isSymmetric() const
Test for symmetry.
void symmetricEigenvaluesEigenvectors(std::vector< T > &eigvals, RankTwoTensorTempl< T > &eigvecs) const
computes eigenvalues and eigenvectors, assuming tens is symmetric, and places them in ascending order...
RankTwoTensorTempl< T > operator-() const
Return the negation of this tensor.
RankFourTensorTempl< T > times(const RankTwoTensorTempl< T > &b) const
Return the general tensor product of this second order tensor and another second order tensor defined...
void syev(const char *calculation_type, std::vector< T > &eigvals, std::vector< T > &a) const
Uses the petscblaslapack.h LAPACKsyev_ routine to find, for symmetric _coords: (1) the eigenvalues (i...
void printDualReal(unsigned int nDual, std::ostream &stm=Moose::out) const
Print the Real part of the RankTwoTensorTempl<ADReal> along with its first nDual dual numbers...
static RankTwoTensorTempl< T > outerProduct(const libMesh::TypeVector< T > &, const libMesh::TypeVector< T > &)
Initialize a second order tensor as the outer product of two vectors, i.e.
IntRange< T > make_range(T beg, T end)
void getRUDecompositionRotation(RankTwoTensorTempl< T > &rot) const
Uses the petscblaslapack.h LAPACKsyev_ routine to perform RU decomposition and obtain the rotation te...
RankTwoTensorTempl< T > ddet() const
Denote the _coords[i][j] by A_ij, then this returns d(det)/dA_ij.
FillMethod
To fill up the 9 entries in the 2nd-order tensor, fillFromInputVector is called with one of the follo...
static RankTwoTensorTempl initializeFromColumns(const libMesh::TypeVector< T > &col0, const libMesh::TypeVector< T > &col1, const libMesh::TypeVector< T > &col2)
Named constructor for initializing from row vectors.
RankTwoTensorTempl< T > initialContraction(const RankFourTensorTempl< T > &b) const
returns this_ij * b_ijkl
void mooseSetToZero< RankTwoTensor >(RankTwoTensor &v)
Helper function template specialization to set an object to zero.
RankTwoTensorTempl(const libMesh::TensorValue< T > &a)
The conversion operator from a libMesh::TensorValue
InitMethod
The initialization method.
RankTwoTensorTempl< typename CompareTypes< T, T2 >::supertype > operator/(const T2 &a) const
Return this tensor divided by a scalar (component-wise)
auto min(const L &left, const R &right)
RankTwoTensorTempl(const SymmetricRankTwoTensorTempl< T2 > &a)
The conversion operator from a SymmetricRankTwoTensorTempl
void fillFromScalarVariable(const VariableValue &scalar_variable)
The smart mutator that determines how to fill the second order tensor based on the order of the scal...
MooseUnits pow(const MooseUnits &, int)
static constexpr Real sqrt2
std::sqrt is not constexpr, so we add sqrt(2) as a constant (used in Mandel notation) ...
MooseArray< Real > VariableValue