Reference documentation for deal.II version GIT relicensing-214-g6e74dec06b 2024-03-27 18:10:01+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
Classes | Namespaces | Enumerations | Functions
symmetric_tensor.h File Reference
#include <deal.II/base/config.h>
#include <deal.II/base/ndarray.h>
#include <deal.II/base/numbers.h>
#include <deal.II/base/table_indices.h>
#include <deal.II/base/template_constraints.h>
#include <deal.II/base/tensor.h>
#include <array>

Go to the source code of this file.

Classes

struct  internal::ProductTypeImpl< SymmetricTensor< rank, dim, T >, std::complex< U > >
 
struct  internal::ProductTypeImpl< SymmetricTensor< rank, dim, std::complex< T > >, std::complex< U > >
 
struct  internal::ProductTypeImpl< std::complex< T >, SymmetricTensor< rank, dim, U > >
 
struct  internal::ProductTypeImpl< std::complex< T >, SymmetricTensor< rank, dim, std::complex< U > > >
 
struct  internal::SymmetricTensorAccessors::double_contraction_result< rank1, rank2, dim, Number, OtherNumber >
 
struct  internal::SymmetricTensorAccessors::double_contraction_result< 2, 2, dim, Number, OtherNumber >
 
struct  internal::SymmetricTensorAccessors::StorageType< 2, dim, Number >
 
struct  internal::SymmetricTensorAccessors::StorageType< 4, dim, Number >
 
struct  internal::SymmetricTensorAccessors::AccessorTypes< rank, dim, true, Number >
 
struct  internal::SymmetricTensorAccessors::AccessorTypes< rank, dim, false, Number >
 
class  internal::SymmetricTensorAccessors::Accessor< rank, dim, constness, P, Number >
 
class  internal::SymmetricTensorAccessors::Accessor< rank, dim, constness, 1, Number >
 
class  SymmetricTensor< rank_, dim, Number >
 
struct  internal::SymmetricTensorImplementation::SortEigenValuesVectors< dim, Number >
 

Namespaces

namespace  internal
 
namespace  internal::SymmetricTensorImplementation
 
namespace  internal::SymmetricTensorAccessors
 

Enumerations

enum struct  SymmetricTensorEigenvectorMethod { hybrid , ql_implicit_shifts , jacobi }
 

Functions

template<int dim, typename Number = double>
DEAL_II_HOST constexpr SymmetricTensor< 2, dim, Number > unit_symmetric_tensor ()
 
template<int dim, typename Number = double>
DEAL_II_HOST constexpr SymmetricTensor< 4, dim, Number > deviator_tensor ()
 
template<int dim, typename Number = double>
DEAL_II_HOST constexpr SymmetricTensor< 4, dim, Number > identity_tensor ()
 
template<int dim, typename Number >
DEAL_II_HOST constexpr SymmetricTensor< 2, dim, Number > invert (const SymmetricTensor< 2, dim, Number > &)
 
template<int dim, typename Number >
DEAL_II_HOST constexpr SymmetricTensor< 4, dim, Number > invert (const SymmetricTensor< 4, dim, Number > &)
 
template<int dim2, typename Number >
DEAL_II_HOST constexpr Number trace (const SymmetricTensor< 2, dim2, Number > &)
 
template<int dim, typename Number >
DEAL_II_HOST constexpr SymmetricTensor< 2, dim, Number > deviator (const SymmetricTensor< 2, dim, Number > &)
 
template<int dim, typename Number >
DEAL_II_HOST constexpr Number determinant (const SymmetricTensor< 2, dim, Number > &)
 
DEAL_II_HOST constexpr TableIndices< 2 > internal::SymmetricTensorAccessors::merge (const TableIndices< 2 > &previous_indices, const unsigned int new_index, const unsigned int position)
 
DEAL_II_HOST constexpr TableIndices< 4 > internal::SymmetricTensorAccessors::merge (const TableIndices< 4 > &previous_indices, const unsigned int new_index, const unsigned int position)
 
template<int rank_, int dim, typename Number , typename OtherNumber >
DEAL_II_HOST constexpr SymmetricTensor< rank_, dim, typename ProductType< Number, OtherNumber >::type > operator+ (const SymmetricTensor< rank_, dim, Number > &left, const SymmetricTensor< rank_, dim, OtherNumber > &right)
 
template<int rank_, int dim, typename Number , typename OtherNumber >
DEAL_II_HOST constexpr SymmetricTensor< rank_, dim, typename ProductType< Number, OtherNumber >::type > operator- (const SymmetricTensor< rank_, dim, Number > &left, const SymmetricTensor< rank_, dim, OtherNumber > &right)
 
template<int rank_, int dim, typename Number , typename OtherNumber >
DEAL_II_HOST constexpr Tensor< rank_, dim, typename ProductType< Number, OtherNumber >::type > operator+ (const SymmetricTensor< rank_, dim, Number > &left, const Tensor< rank_, dim, OtherNumber > &right)
 
template<int rank_, int dim, typename Number , typename OtherNumber >
DEAL_II_HOST constexpr Tensor< rank_, dim, typename ProductType< Number, OtherNumber >::type > operator+ (const Tensor< rank_, dim, Number > &left, const SymmetricTensor< rank_, dim, OtherNumber > &right)
 
template<int rank_, int dim, typename Number , typename OtherNumber >
DEAL_II_HOST constexpr Tensor< rank_, dim, typename ProductType< Number, OtherNumber >::type > operator- (const SymmetricTensor< rank_, dim, Number > &left, const Tensor< rank_, dim, OtherNumber > &right)
 
template<int rank_, int dim, typename Number , typename OtherNumber >
DEAL_II_HOST constexpr Tensor< rank_, dim, typename ProductType< Number, OtherNumber >::type > operator- (const Tensor< rank_, dim, Number > &left, const SymmetricTensor< rank_, dim, OtherNumber > &right)
 
template<int dim, typename Number >
DEAL_II_HOST constexpr Number third_invariant (const SymmetricTensor< 2, dim, Number > &t)
 
template<int dim, typename Number >
DEAL_II_HOST constexpr Number trace (const SymmetricTensor< 2, dim, Number > &d)
 
template<int dim, typename Number >
DEAL_II_HOST constexpr Number first_invariant (const SymmetricTensor< 2, dim, Number > &t)
 
template<typename Number >
DEAL_II_HOST constexpr Number second_invariant (const SymmetricTensor< 2, 1, Number > &)
 
template<typename Number >
DEAL_II_HOST constexpr Number second_invariant (const SymmetricTensor< 2, 2, Number > &t)
 
template<typename Number >
DEAL_II_HOST constexpr Number second_invariant (const SymmetricTensor< 2, 3, Number > &t)
 
template<typename Number >
std::array< Number, 1 > eigenvalues (const SymmetricTensor< 2, 1, Number > &T)
 
template<typename Number >
std::array< Number, 2 > eigenvalues (const SymmetricTensor< 2, 2, Number > &T)
 
template<typename Number >
std::array< Number, 3 > eigenvalues (const SymmetricTensor< 2, 3, Number > &T)
 
template<int dim, typename Number >
void internal::SymmetricTensorImplementation::tridiagonalize (const ::SymmetricTensor< 2, dim, Number > &A, ::Tensor< 2, dim, Number > &Q, std::array< Number, dim > &d, std::array< Number, dim - 1 > &e)
 
template<int dim, typename Number >
std::array< std::pair< Number, Tensor< 1, dim, Number > >, dim > internal::SymmetricTensorImplementation::ql_implicit_shifts (const ::SymmetricTensor< 2, dim, Number > &A)
 
template<int dim, typename Number >
std::array< std::pair< Number, Tensor< 1, dim, Number > >, dim > internal::SymmetricTensorImplementation::jacobi (::SymmetricTensor< 2, dim, Number > A)
 
template<typename Number >
std::array< std::pair< Number, Tensor< 1, 2, Number > >, 2 > internal::SymmetricTensorImplementation::hybrid (const ::SymmetricTensor< 2, 2, Number > &A)
 
template<typename Number >
std::array< std::pair< Number, Tensor< 1, 3, Number > >, 3 > internal::SymmetricTensorImplementation::hybrid (const ::SymmetricTensor< 2, 3, Number > &A)
 
template<int dim, typename Number >
std::array< std::pair< Number, Tensor< 1, dim, Number > >, std::integral_constant< int, dim >::value > eigenvectors (const SymmetricTensor< 2, dim, Number > &T, const SymmetricTensorEigenvectorMethod method=SymmetricTensorEigenvectorMethod::ql_implicit_shifts)
 
template<int rank_, int dim, typename Number >
DEAL_II_HOST constexpr SymmetricTensor< rank_, dim, Number > transpose (const SymmetricTensor< rank_, dim, Number > &t)
 
template<int dim, typename Number >
DEAL_II_HOST constexpr SymmetricTensor< 4, dim, Number > outer_product (const SymmetricTensor< 2, dim, Number > &t1, const SymmetricTensor< 2, dim, Number > &t2)
 
template<int dim, typename Number >
std::pair< SymmetricTensor< 2, dim, Number >, SymmetricTensor< 2, dim, Number > > positive_negative_split (const SymmetricTensor< 2, dim, Number > &original_tensor)
 
template<int dim, typename Number >
std::tuple< SymmetricTensor< 2, dim, Number >, SymmetricTensor< 2, dim, Number >, SymmetricTensor< 4, dim, Number >, SymmetricTensor< 4, dim, Number > > positive_negative_projectors (const SymmetricTensor< 2, dim, Number > &original_tensor)
 
template<int dim, typename Number >
DEAL_II_HOST constexpr SymmetricTensor< 2, dim, Number > symmetrize (const Tensor< 2, dim, Number > &t)
 
template<int rank_, int dim, typename Number >
DEAL_II_HOST constexpr SymmetricTensor< rank_, dim, Number > operator* (const SymmetricTensor< rank_, dim, Number > &t, const Number &factor)
 
template<int rank_, int dim, typename Number >
DEAL_II_HOST constexpr SymmetricTensor< rank_, dim, Number > operator* (const Number &factor, const SymmetricTensor< rank_, dim, Number > &t)
 

Enumeration Type Documentation

◆ SymmetricTensorEigenvectorMethod

An enumeration for the algorithm to be employed when performing the computation of normalized eigenvectors and their corresponding eigenvalues by the eigenvalues() and eigenvectors() methods operating on SymmetricTensor objects.

The specialized algorithms utilized in computing the eigenvectors are presented in

@article{Kopp2008,
title = {Efficient numerical diagonalization of hermitian 3x3
matrices},
author = {Kopp, J.},
journal = {International Journal of Modern Physics C},
year = {2008},
volume = {19},
number = {3},
pages = {523--548},
doi = {10.1142/S0129183108012303},
eprinttype = {arXiv},
eprint = {physics/0610206v3},
eprintclass = {physics.comp-ph},
url =
{https://www.mpi-hd.mpg.de/personalhomes/globes/3x3/index.html}
}
Enumerator
hybrid 

A hybrid approach that preferentially uses the characteristic equation to compute eigenvalues and an analytical approach based on the cross-product for the eigenvectors. If the computations are deemed too inaccurate then the method falls back to ql_implicit_shifts.

This method potentially offers the quickest computation if the pathological case is not encountered.

ql_implicit_shifts 

The iterative QL algorithm with implicit shifts applied after tridiagonalization of the tensor using the householder method.

This method offers a compromise between speed of computation and its robustness. This method is particularly useful when the elements of \(T\) have greatly varying magnitudes, which would typically lead to a loss of accuracy when computing the smaller eigenvalues.

jacobi 

The iterative Jacobi algorithm.

This method offers is the most robust of the available options, with reliable results obtained for even the most pathological cases. It is, however, the slowest algorithm of all of those implemented.

Definition at line 3182 of file symmetric_tensor.h.

Function Documentation

◆ unit_symmetric_tensor()

template<int dim, typename Number = double>
DEAL_II_HOST constexpr SymmetricTensor< 2, dim, Number > unit_symmetric_tensor ( )
inlineconstexpr

Return a unit symmetric tensor of rank 2, i.e., the \(\text{dim}\times\text{dim}\) identity matrix \(\mathbf I\). For example, if dim==2, then this matrix has the form

\[ I_{2\times 2} = \begin{pmatrix} 1 & 0 \\ 0 & 1 \end{pmatrix}. \]

Definition at line 3290 of file symmetric_tensor.h.

◆ deviator_tensor()

template<int dim, typename Number = double>
DEAL_II_HOST constexpr SymmetricTensor< 4, dim, Number > deviator_tensor ( )
inlineconstexpr

Return the tensor of rank 4 that, when multiplied by a symmetric rank 2 tensor \(\mathbf T\) returns the deviator \(\text{dev}\ \mathbf T\). It is the operator representation of the linear deviator operator \(\mathbb P\), also known as the volumetric projection tensor, calculated as:

\begin{align*} \mathbb{P} &=\mathbb{S} -\frac{1}{\text{dim}} \mathbf I \otimes \mathbf I \\ \mathcal{P}_{ijkl} &= \frac 12 \left(\delta_{ik} \delta_{jl} + \delta_{il} \delta_{jk} \right) - \frac{1}{\text{dim}} \delta_{ij} \delta_{kl} \end{align*}

For every tensor T, there holds the identity deviator<dim,Number>(T) == deviator_tensor<dim,Number>() * T, up to numerical round-off.

\[ \text{dev}\mathbf T = \mathbb P : \mathbf T \]

Note
The reason this operator representation is provided is to simplify taking derivatives of the deviatoric part of tensors:

\[ \frac{\partial \text{dev}\mathbf{T}}{\partial \mathbf T} = \mathbb P. \]

Definition at line 3318 of file symmetric_tensor.h.

◆ identity_tensor()

template<int dim, typename Number = double>
DEAL_II_HOST constexpr SymmetricTensor< 4, dim, Number > identity_tensor ( )
inlineconstexpr

Return the fourth-order symmetric identity tensor \(\mathbb S\) which maps symmetric second-order tensors, such as \(\mathbf A\), to themselves.

\[ \mathbb S : \mathbf A = \mathbf A \]

Note that this tensor, even though it is the identity, has a somewhat funny form, and in particular does not only consist of zeros and ones. For example, for dim=2, the identity tensor has all zero entries except for

\[ \mathcal{S}_{0000} = \mathcal{S}_{1111} = 1 \]

\[ \mathcal{S}_{0101} = \mathcal{S}_{0110} = \mathcal{S}_{1001} = \mathcal{S}_{1010} = \frac 12. \]

In index notation, we can write the general form

\[ \mathcal{S}_{ijkl} = \frac 12 \left( \delta_{ik} \delta_{jl} + \delta_{il} \delta_{jk} \right). \]

To see why this factor of \(1 / 2\) is necessary, consider computing \(\mathbf A= \mathbb I : \mathbf B\). For the element \(A_{01}\) we have \(A_{01} = \mathcal{S}_{0100} B_{00} + \mathcal{S}_{0111} B_{11} + \mathcal{S}_{0101} B_{01} + \mathcal{S}_{0110} B_{10}\). On the other hand, we need to have \(A_{01} = B_{01}\), and symmetry implies \(B_{01}=B_{10}\), leading to \(A_{01} = (\mathcal{S}_{0101} + \mathcal{S}_{0110}) B_{01}\), or, again by symmetry, \(\mathcal{S}_{0101} = \mathcal{S}_{0110} = \frac 12\). Similar considerations hold for the three-dimensional case.

This issue is also explained in the introduction to step-44.

Definition at line 3347 of file symmetric_tensor.h.

◆ invert() [1/2]

template<int dim, typename Number >
DEAL_II_HOST constexpr SymmetricTensor< 2, dim, Number > invert ( const SymmetricTensor< 2, dim, Number > &  t)
constexpr

Invert a symmetric rank-2 tensor.

Note
If a tensor is not invertible, then the result is unspecified, but will likely contain the results of a division by zero or a very small number at the very least.

Definition at line 3382 of file symmetric_tensor.h.

◆ invert() [2/2]

template<int dim, typename Number >
DEAL_II_HOST constexpr SymmetricTensor< 4, dim, Number > invert ( const SymmetricTensor< 4, dim, Number > &  t)
constexpr

Invert a symmetric rank-4 tensor. Since symmetric rank-4 tensors are mappings from and to symmetric rank-2 tensors, they can have an inverse.

If a tensor is not invertible, then the result is unspecified, but will likely contain the results of a division by zero or a very small number at the very least.

Definition at line 3402 of file symmetric_tensor.h.

◆ trace() [1/2]

template<int dim2, typename Number >
DEAL_II_HOST constexpr Number trace ( const SymmetricTensor< 2, dim2, Number > &  )
inlineconstexpr

Compute and return the trace of a tensor of rank 2, i.e. the sum of its diagonal entries. The trace is the first invariant of a rank-2 tensor.

\[ \text{tr} \mathbf A = \sum_i A_{ii} \]

◆ deviator()

template<int dim, typename Number >
DEAL_II_HOST constexpr SymmetricTensor< 2, dim, Number > deviator ( const SymmetricTensor< 2, dim, Number > &  t)
inlineconstexpr

Compute the deviator of a symmetric tensor, which is defined as \(\text{dev} \mathbf T = \mathbf T - \frac{1}{\text{dim}} \text{tr}\mathbf T \; \mathbf I\), where \(\mathbf I\) is the identity operator. This quantity equals the original tensor minus its contractive or dilative component and refers to the shear in, for example, elasticity.

Definition at line 3273 of file symmetric_tensor.h.

◆ determinant()

template<int dim, typename Number >
DEAL_II_HOST constexpr Number determinant ( const SymmetricTensor< 2, dim, Number > &  t)
inlineconstexpr

Compute the determinant of a rank 2 symmetric tensor. The determinant is also commonly referred to as the third invariant of rank-2 tensors.

For a one-dimensional tensor, the determinant equals the only element and is therefore equivalent to the trace.

For greater notational simplicity, there is also a third_invariant() function that returns the determinant of a tensor.

Definition at line 2724 of file symmetric_tensor.h.

◆ operator+() [1/3]

template<int rank_, int dim, typename Number , typename OtherNumber >
DEAL_II_HOST constexpr SymmetricTensor< rank_, dim, typename ProductType< Number, OtherNumber >::type > operator+ ( const SymmetricTensor< rank_, dim, Number > &  left,
const SymmetricTensor< rank_, dim, OtherNumber > &  right 
)
inlineconstexpr

Addition of two symmetric tensors of equal rank. The result is another SymmetricTensor that has a number type that is compatible with the operation.

If possible (e.g. when Number and OtherNumber are of the same type, or if the result of Number() + OtherNumber() is another Number), you should use operator+= instead since this does not require the creation of a temporary variable.

Definition at line 2618 of file symmetric_tensor.h.

◆ operator-() [1/3]

template<int rank_, int dim, typename Number , typename OtherNumber >
DEAL_II_HOST constexpr SymmetricTensor< rank_, dim, typename ProductType< Number, OtherNumber >::type > operator- ( const SymmetricTensor< rank_, dim, Number > &  left,
const SymmetricTensor< rank_, dim, OtherNumber > &  right 
)
inlineconstexpr

Subtraction of two symmetric tensors of equal rank. The result is another SymmetricTensor that has a number type that is compatible with the operation.

If possible (e.g. when Number and OtherNumber are of the same type, or if the result of Number() - OtherNumber() is another Number), you should use operator-= instead since this does not require the creation of a temporary variable.

Definition at line 2643 of file symmetric_tensor.h.

◆ operator+() [2/3]

template<int rank_, int dim, typename Number , typename OtherNumber >
DEAL_II_HOST constexpr Tensor< rank_, dim, typename ProductType< Number, OtherNumber >::type > operator+ ( const SymmetricTensor< rank_, dim, Number > &  left,
const Tensor< rank_, dim, OtherNumber > &  right 
)
constexpr

Addition of a SymmetricTensor and a general Tensor of equal rank. The result is a general Tensor that has a number type that is compatible with the operation.

Definition at line 2663 of file symmetric_tensor.h.

◆ operator+() [3/3]

template<int rank_, int dim, typename Number , typename OtherNumber >
DEAL_II_HOST constexpr Tensor< rank_, dim, typename ProductType< Number, OtherNumber >::type > operator+ ( const Tensor< rank_, dim, Number > &  left,
const SymmetricTensor< rank_, dim, OtherNumber > &  right 
)
constexpr

Addition of a general Tensor with a SymmetricTensor of equal rank. The result is a general Tensor that has a number type that is compatible with the operation.

Definition at line 2680 of file symmetric_tensor.h.

◆ operator-() [2/3]

template<int rank_, int dim, typename Number , typename OtherNumber >
DEAL_II_HOST constexpr Tensor< rank_, dim, typename ProductType< Number, OtherNumber >::type > operator- ( const SymmetricTensor< rank_, dim, Number > &  left,
const Tensor< rank_, dim, OtherNumber > &  right 
)
constexpr

Subtraction of a general Tensor from a SymmetricTensor of equal rank. The result is a general Tensor that has a number type that is compatible with the operation.

Definition at line 2697 of file symmetric_tensor.h.

◆ operator-() [3/3]

template<int rank_, int dim, typename Number , typename OtherNumber >
DEAL_II_HOST constexpr Tensor< rank_, dim, typename ProductType< Number, OtherNumber >::type > operator- ( const Tensor< rank_, dim, Number > &  left,
const SymmetricTensor< rank_, dim, OtherNumber > &  right 
)
constexpr

Subtraction of a SymmetricTensor from a general Tensor of equal rank. The result is a general Tensor that has a number type that is compatible with the operation.

Definition at line 2714 of file symmetric_tensor.h.

◆ third_invariant()

template<int dim, typename Number >
DEAL_II_HOST constexpr Number third_invariant ( const SymmetricTensor< 2, dim, Number > &  t)
constexpr

Compute the determinant of a rank 2 symmetric tensor. This function therefore computes the same value as the determinant() functions and is only provided for greater notational simplicity (since there are also functions first_invariant() and second_invariant()).

\[ I_3 (\mathbf A) = III (\mathbf A) = \det (\mathbf A) \]

Definition at line 2764 of file symmetric_tensor.h.

◆ trace() [2/2]

template<int dim, typename Number >
DEAL_II_HOST constexpr Number trace ( const SymmetricTensor< 2, dim, Number > &  d)
inlineconstexpr

Definition at line 2773 of file symmetric_tensor.h.

◆ first_invariant()

template<int dim, typename Number >
DEAL_II_HOST constexpr Number first_invariant ( const SymmetricTensor< 2, dim, Number > &  t)
constexpr

Compute the trace of a rank 2 symmetric tensor. This function therefore computes the same value as the trace() functions and is only provided for greater notational simplicity (since there are also functions second_invariant() and third_invariant()).

\[ I_1 (\mathbf A) = I (\mathbf A) = \text{tr} \mathbf A = \sum_i A_{ii} \]

Definition at line 2795 of file symmetric_tensor.h.

◆ second_invariant() [1/3]

template<typename Number >
DEAL_II_HOST constexpr Number second_invariant ( const SymmetricTensor< 2, 1, Number > &  )
constexpr

Compute the second invariant of a tensor of rank 2. The second invariant of a tensor \(\mathbf A\) is defined as \(I_2 (\mathbf A) = II(\mathbf A) = \frac 12 \left[ (\text{tr} \mathbf A)^2 - \text{tr} (\mathbf{A}^2) \right]\).

For the kind of arguments to this function, i.e., a rank-2 tensor of size 1, the result is simply zero.

Definition at line 2814 of file symmetric_tensor.h.

◆ second_invariant() [2/3]

template<typename Number >
DEAL_II_HOST constexpr Number second_invariant ( const SymmetricTensor< 2, 2, Number > &  t)
constexpr

Compute the second invariant of a tensor of rank 2. The second invariant of a tensor \(\mathbf A\) is defined as \(I_2 (\mathbf A) = II(\mathbf A) = \frac 12 \left[ (\text{tr} \mathbf A)^2 - \text{tr} (\mathbf{A}^2) \right]\).

For the kind of arguments to this function, i.e., a symmetric rank-2 tensor of size 2, the result is (counting indices starting at one) \(I_2(\mathbf A) = II(\mathbf A) = \frac 12 \left[ (A_{11} + A_{22})^2 - (A_{11}^2+2 A_{12}^2+ A_{22}^2) \right] = A_{11} A_{22} - A_{12}^2\). As expected, for the \(2\times 2\) symmetric tensors this function handles, this equals the determinant of the tensor. (This is so because for \(2\times 2\) symmetric tensors, there really are only two invariants, so the second and third invariant are the same; the determinant is the third invariant.)

Definition at line 2841 of file symmetric_tensor.h.

◆ second_invariant() [3/3]

template<typename Number >
DEAL_II_HOST constexpr Number second_invariant ( const SymmetricTensor< 2, 3, Number > &  t)
constexpr

Compute the second invariant of a tensor of rank 2. The second invariant of a tensor \(\mathbf A\) is defined as \(I_2 (\mathbf A) = II(\mathbf A) = \frac 12 \left[ (\text{tr} \mathbf A)^2 - \text{tr} (\mathbf{A}^2) \right]\).

Definition at line 2858 of file symmetric_tensor.h.

◆ eigenvalues() [1/3]

template<typename Number >
std::array< Number, 1 > eigenvalues ( const SymmetricTensor< 2, 1, Number > &  T)

Return the eigenvalues of a symmetric \(1 \times 1\) tensor. The (single) entry of the tensor is, of course, equal to the (single) eigenvalue.

◆ eigenvalues() [2/3]

template<typename Number >
std::array< Number, 2 > eigenvalues ( const SymmetricTensor< 2, 2, Number > &  T)

Return the eigenvalues of a symmetric \(2\times 2\) tensor. The array of eigenvalues is sorted in descending order.

For \(2\times 2\) tensors, the eigenvalues of tensor \(\mathbf T\) are the roots of the characteristic polynomial \(0 = \lambda^2 - \lambda\;\text{tr}\mathbf{T} + \det \mathbf{T}\) as given by \(\lambda_1, \lambda_2 = \frac{1}{2} \left[ \text{tr} \mathbf{T} \pm \sqrt{(\text{tr} \mathbf{T})^2 - 4 \det \mathbf{T}} \right]\).

Warning
The algorithm employed here determines the eigenvalues by computing the roots of the characteristic polynomial. In the case that there exists a common root (the eigenvalues are equal), the computation is subject to round-off errors of order \(\sqrt{\epsilon}\). As an alternative, the eigenvectors() function provides a more robust, but costly, method to compute the eigenvalues of a symmetric tensor.

◆ eigenvalues() [3/3]

template<typename Number >
std::array< Number, 3 > eigenvalues ( const SymmetricTensor< 2, 3, Number > &  T)

Return the eigenvalues of a symmetric \(3\times 3\) tensor. The array of eigenvalues is sorted in descending order.

For \(3\times 3\) tensors, the eigenvalues of tensor \(\mathbf T\) are the roots of the characteristic polynomial \(0 = \lambda^3 - \lambda^2\;\text{tr}\mathbf T - \frac{1}{2} \lambda \left[\text{tr}(\mathbf{T}^2) - (\text{tr}\mathbf T)^2\right] - \det \mathbf T\).

Warning
The algorithm employed here determines the eigenvalues by computing the roots of the characteristic polynomial. In the case that there exists a common root (the eigenvalues are equal), the computation is subject to round-off errors of order \(\sqrt{\epsilon}\). As an alternative, the eigenvectors() function provides a more robust, but costly, method to compute the eigenvalues of a symmetric tensor.

◆ eigenvectors()

template<int dim, typename Number >
std::array< std::pair< Number, Tensor< 1, dim, Number > >, std::integral_constant< int, dim >::value > eigenvectors ( const SymmetricTensor< 2, dim, Number > &  T,
const SymmetricTensorEigenvectorMethod  method = SymmetricTensorEigenvectorMethod::ql_implicit_shifts 
)

Return the eigenvalues and eigenvectors of a real-valued rank-2 symmetric tensor \(\mathbf T\). The array of matched eigenvalue and eigenvector pairs is sorted in descending order (determined by the eigenvalues).

The specialized algorithms utilized in computing the eigenvectors are presented in

@article{Kopp2008,
title = {Efficient numerical diagonalization of hermitian 3x3
matrices},
author = {Kopp, J.},
journal = {International Journal of Modern Physics C},
year = {2008},
volume = {19},
number = {3},
pages = {523--548},
doi = {10.1142/S0129183108012303},
eprinttype = {arXiv},
eprint = {physics/0610206v3},
eprintclass = {physics.comp-ph},
url =
{https://www.mpi-hd.mpg.de/personalhomes/globes/3x3/index.html}
}

◆ transpose()

template<int rank_, int dim, typename Number >
DEAL_II_HOST constexpr SymmetricTensor< rank_, dim, Number > transpose ( const SymmetricTensor< rank_, dim, Number > &  t)
constexpr

Return the transpose of the given symmetric tensor. Since we are working with symmetric objects, the transpose is of course the same as the original tensor. This function mainly exists for compatibility with the Tensor class.

Definition at line 3263 of file symmetric_tensor.h.

◆ outer_product()

template<int dim, typename Number >
DEAL_II_HOST constexpr SymmetricTensor< 4, dim, Number > outer_product ( const SymmetricTensor< 2, dim, Number > &  t1,
const SymmetricTensor< 2, dim, Number > &  t2 
)
inlineconstexpr

Return the tensor of rank 4 that is the outer product of the two tensors given as arguments, i.e. the result \(\mathbb A = \mathbf{T}_1 \otimes \mathbf{T}_2\) satisfies \(\mathbb A : \mathbf B = (\mathbf{T}_2 : \mathbf B) \mathbf{T}_1\) for all symmetric tensors \(\mathbf B\). In index notation

\[ \mathcal{A}_{ijkl} = (T_1)_{ij} (T_2)_{kl} \]

For example, the deviator tensor \(\mathbb P = \mathbb I - \frac{1}{\text{dim}} \mathbf I \otimes \mathbf I\) can be computed as identity_tensor<dim>() - 1/d * outer_product (unit_symmetric_tensor<dim>(), unit_symmetric_tensor<dim>()), since the (double) contraction with the unit tensor yields the trace of a symmetric tensor ( \(\mathbf I : \mathbf B = \text{tr} \mathbf B\)).

Definition at line 3433 of file symmetric_tensor.h.

◆ positive_negative_split()

template<int dim, typename Number >
std::pair< SymmetricTensor< 2, dim, Number >, SymmetricTensor< 2, dim, Number > > positive_negative_split ( const SymmetricTensor< 2, dim, Number > &  original_tensor)

Perform a spectrum decomposition of a 2nd-order symmetric tensor original_tensor given as the input argument,

\[ \mathrm{original\_tensor} = \sum_i \lambda_i \, \boldsymbol{n}_i \otimes \boldsymbol{n}_i, \]

where \(\lambda_i\) is the eigenvalue, and \(\boldsymbol{n}_i\) is the corresponding eigenvector. The output is a pair of 2nd-order symmetric tensors. The first term in the pair is the positive part of the input tensor, and the second term in the pair is the negative part of the input tensor, that is,

\[ \mathrm{positive\_part\_tensor} = \sum_i \left<\lambda_i\right>_+ \boldsymbol{n}_i \otimes \boldsymbol{n}_i, \quad \mathrm{negative\_part\_tensor} = \sum_i \left<\lambda_i\right>_- \boldsymbol{n}_i \otimes \boldsymbol{n}_i, \]

where \(\left<\lambda_i\right>_+ = \mathrm{max}\{ \lambda_i, 0 \}\) and \(\left<\lambda_i\right>_- = \mathrm{min}\{ \lambda_i, 0 \}\). Obviously,

\[ \mathrm{positive\_part\_tensor} + \mathrm{negative\_part\_tensor} = \mathrm{original\_tensor}. \]

Parameters
[in]original_tensorThe 2nd-order symmetric tensor to be split into the positive and negative parts

Definition at line 3475 of file symmetric_tensor.h.

◆ positive_negative_projectors()

template<int dim, typename Number >
std::tuple< SymmetricTensor< 2, dim, Number >, SymmetricTensor< 2, dim, Number >, SymmetricTensor< 4, dim, Number >, SymmetricTensor< 4, dim, Number > > positive_negative_projectors ( const SymmetricTensor< 2, dim, Number > &  original_tensor)

This function is similar to the function positive_negative_split(). That is, perform a spectrum decomposition of a 2nd-order symmetric tensor original_tensor given as the input argument, and split it into a positive part and a negative part. Moreover, this function also provides the derivatives. Let \(\mathbf{A}\) represent the input 2nd-order symmetric tensor original_tensor, \(\mathbf{A}^+\) represent the positive part, and \(\mathbf{A}^-\) represent the negative part. Then, two fourth-order tensors are defined as

\[ \mathbb{P}^+ = \frac{\partial \mathbf{A}^+}{\partial \mathbf{A}}, \quad \mathbb{P}^- = \frac{\partial \mathbf{A}^-}{\partial \mathbf{A}}, \]

where \(\mathbb{P}^+\) is the positive projector and \(\mathbb{P}^-\) is the negative projector. These two fourth-order tensors satisfy the following properties:

\[ \mathbb{P}^+ : \mathbf{A} = \mathbf{A}^+, \quad \mathbb{P}^- : \mathbf{A} = \mathbf{A}^-. \]

Since \(\mathbb{P}^+\) and \(\mathbb{P}^-\) are 4th-order projectors,

\[ \mathbb{P}^+ : \mathbf{A}^+ = \mathbf{A}^+, \quad \mathbb{P}^- : \mathbf{A}^- = \mathbf{A}^-, \quad \mathbb{P}^+ : \mathbf{A}^- = \mathbb{P}^- : \mathbf{A}^+ = \mathbf{0}. \]

Lastly,

\[ \mathbb{P}^+ + \mathbb{P}^- = \mathbb{S}, \]

where \(\mathbb{S}\) is the fourth-order symmetric identity tensor Physics::Elasticity::StandardTensors< dim >::S. The output of this function is a tuple containing four terms. The first term is \(\mathbf{A}^+\), the second term is \(\mathbf{A}^-\), the third term is \(\mathbb{P}^+\), and the fourth term is \(\mathbb{P}^-\).

Parameters
[in]original_tensorThe 2nd-order symmetric tensor to be split into the positive and negative parts

Definition at line 3541 of file symmetric_tensor.h.

◆ symmetrize()

template<int dim, typename Number >
DEAL_II_HOST constexpr SymmetricTensor< 2, dim, Number > symmetrize ( const Tensor< 2, dim, Number > &  t)
inlineconstexpr

Return the symmetrized version of a full rank-2 tensor, i.e. \(\text{sym}\mathbf A = \frac 12 \left(\mathbf A + \mathbf{A}^T\right)\), as a symmetric rank-2 tensor. This is the version for general dimensions.

Definition at line 3665 of file symmetric_tensor.h.

◆ operator*() [1/2]

template<int rank_, int dim, typename Number >
DEAL_II_HOST constexpr SymmetricTensor< rank_, dim, Number > operator* ( const SymmetricTensor< rank_, dim, Number > &  t,
const Number &  factor 
)
inlineconstexpr

Multiplication of a symmetric tensor of general rank with a scalar from the right. This version of the operator is used if the scalar has the same data type as is used to store the elements of the symmetric tensor.

Definition at line 3690 of file symmetric_tensor.h.

◆ operator*() [2/2]

template<int rank_, int dim, typename Number >
DEAL_II_HOST constexpr SymmetricTensor< rank_, dim, Number > operator* ( const Number &  factor,
const SymmetricTensor< rank_, dim, Number > &  t 
)
constexpr

Multiplication of a symmetric tensor of general rank with a scalar from the left. This version of the operator is used if the scalar has the same data type as is used to store the elements of the symmetric tensor.

Definition at line 3708 of file symmetric_tensor.h.