RankFourTensorTempl is designed to handle any N-dimensional fourth order tensor, C. More...

#include <RankFourTensor.h>

Classes
struct	TwoTensorMultTraits

struct	TwoTensorMultTraits< RankTwoTensorTempl, Scalar >

struct	TwoTensorMultTraits< TensorValue, Scalar >

struct	TwoTensorMultTraits< TypeTensor, Scalar >

Public Types
enum	InitMethod { initNone , initIdentity , initIdentityFour , initIdentitySymmetricFour , initIdentityDeviatoric }
	Initialization method. More...

enum	FillMethod { antisymmetric , symmetric9 , symmetric21 , general_isotropic , symmetric_isotropic , symmetric_isotropic_E_nu , antisymmetric_isotropic , axisymmetric_rz , general , principal , orthotropic }
	To fill up the 81 entries in the 4th-order tensor, fillFromInputVector is called with one of the following fill_methods. More...

typedef tuple_of< 4, unsigned int >	index_type

typedef T	value_type

Public Member Functions
	RankFourTensorTempl ()
	Default constructor; fills to zero.

	RankFourTensorTempl (const InitMethod)
	Select specific initialization pattern.

	RankFourTensorTempl (const std::vector< T > &, FillMethod)
	Fill from vector.

	RankFourTensorTempl (const RankFourTensorTempl< T > &a)=default
	Copy assignment operator must be defined if used.

template<typename T2 >
	RankFourTensorTempl (const RankFourTensorTempl< T2 > &copy)
	Copy constructor.

template<typename T2 >
	RankFourTensorTempl (const SymmetricRankFourTensorTempl< T2 > &t)
	The conversion operator from a `SymmetricRankFourTensorTempl`

T &	operator() (unsigned int i, unsigned int j, unsigned int k, unsigned int l)
	Gets the value for the indices specified. Takes indices ranging from 0-2 for i, j, k, and l.

const T &	operator() (unsigned int i, unsigned int j, unsigned int k, unsigned int l) const
	Gets the value for the indices specified.

void	zero ()
	Zeros out the tensor.

void	print (std::ostream &stm=Moose::out) const
	Print the rank four tensor.

void	printReal (std::ostream &stm=Moose::out) const
	Print the values of the rank four tensor.

RankFourTensorTempl< T > &	operator= (const RankFourTensorTempl< T > &a)=default
	copies values from a into this tensor

template<typename Scalar >
std::enable_if< libMesh::ScalarTraits< Scalar >::value, RankFourTensorTempl & >::type	operator= (const Scalar &libmesh_dbg_var(p))
	Assignment-from-scalar operator.

template<template< typename > class Tensor, typename T2 >
auto	operator* (const Tensor< T2 > &a) const -> typename std::enable_if< TwoTensorMultTraits< Tensor, T2 >::value, RankTwoTensorTempl< decltype(T() *T2())> >::type
	C_ijkl*a_kl.

template<typename T2 >
auto	operator* (const T2 &a) const -> typename std::enable_if< libMesh::ScalarTraits< T2 >::value, RankFourTensorTempl< decltype(T() *T2())> >::type
	C_ijkl*a.

RankFourTensorTempl< T > &	operator*= (const T &a)
	C_ijkl *= a.

template<typename T2 >
auto	operator/ (const T2 &a) const -> typename std::enable_if< libMesh::ScalarTraits< T2 >::value, RankFourTensorTempl< decltype(T()/T2())> >::type
	C_ijkl/a.

RankFourTensorTempl< T > &	operator/= (const T &a)
	C_ijkl /= a for all i, j, k, l.

RankFourTensorTempl< T > &	operator+= (const RankFourTensorTempl< T > &a)
	C_ijkl += a_ijkl for all i, j, k, l.

template<typename T2 >
auto	operator+ (const RankFourTensorTempl< T2 > &a) const -> RankFourTensorTempl< decltype(T()+T2())>
	C_ijkl + a_ijkl.

RankFourTensorTempl< T > &	operator-= (const RankFourTensorTempl< T > &a)
	C_ijkl -= a_ijkl.

template<typename T2 >
auto	operator- (const RankFourTensorTempl< T2 > &a) const -> RankFourTensorTempl< decltype(T() - T2())>
	C_ijkl - a_ijkl.

RankFourTensorTempl< T >	operator- () const
	-C_ijkl

template<typename T2 >
auto	operator* (const RankFourTensorTempl< T2 > &a) const -> RankFourTensorTempl< decltype(T() *T2())>
	C_ijpq*a_pqkl.

T	L2norm () const
	sqrt(C_ijkl*C_ijkl)

RankFourTensorTempl< T >	invSymm () const
	This returns A_ijkl such that C_ijklA_klmn = 0.5(de_im de_jn + de_in de_jm) This routine assumes that C_ijkl = C_jikl = C_ijlk.

RankFourTensorTempl< T >	inverse () const
	This returns A_ijkl such that C_ijkl*A_klmn = de_im de_jn i.e.

void	rotate (const TypeTensor< T > &R)
	Rotate the tensor using C_ijkl = R_im R_jn R_ko R_lp C_mnop.

RankFourTensorTempl< T >	transposeMajor () const
	Transpose the tensor by swapping the first pair with the second pair of indices.

RankFourTensorTempl< T >	transposeIj () const
	Transpose the tensor by swapping the first two indeces.

RankFourTensorTempl< T >	transposeKl () const
	Transpose the tensor by swapping the last two indeces.

template<int m>
RankThreeTensorTempl< T >	contraction (const libMesh::VectorValue< T > &b) const
	single contraction of a RankFourTensor with a vector over index m

void	surfaceFillFromInputVector (const std::vector< T > &input)
	Fills the tensor entries ignoring the last dimension (ie, C_ijkl=0 if any of i, j, k, or l = 3).

void	fillFromInputVector (const std::vector< T > &input, FillMethod fill_method)
	fillFromInputVector takes some number of inputs to fill the Rank-4 tensor.

template<typename T2 >
void	fillSymmetric9FromInputVector (const T2 &input)
	fillSymmetric9FromInputVector takes 9 inputs to fill in the Rank-4 tensor with the appropriate crystal symmetries maintained.

template<typename T2 >
void	fillSymmetric21FromInputVector (const T2 &input)
	fillSymmetric21FromInputVector takes either 21 inputs to fill in the Rank-4 tensor with the appropriate crystal symmetries maintained.

RankTwoTensorTempl< T >	innerProductTranspose (const RankTwoTensorTempl< T > &) const
	Inner product of the major transposed tensor with a rank two tensor.

T	contractionIj (unsigned int, unsigned int, const RankTwoTensorTempl< T > &) const
	Sum C_ijkl M_kl for a given i,j.

T	contractionKl (unsigned int, unsigned int, const RankTwoTensorTempl< T > &) const
	Sum M_ij C_ijkl for a given k,l.

T	sum3x3 () const
	Calculates the sum of Ciijj for i and j varying from 0 to 2.

libMesh::VectorValue< T >	sum3x1 () const
	Calculates the vector a[i] = sum over j Ciijj for i and j varying from 0 to 2.

RankFourTensorTempl< T >	tripleProductJkl (const RankTwoTensorTempl< T > &, const RankTwoTensorTempl< T > &, const RankTwoTensorTempl< T > &) const
	Calculates C_imnt A_jm B_kn C_lt.

RankFourTensorTempl< T >	tripleProductIkl (const RankTwoTensorTempl< T > &, const RankTwoTensorTempl< T > &, const RankTwoTensorTempl< T > &) const
	Calculates C_mjnt A_im B_kn C_lt.

RankFourTensorTempl< T >	tripleProductIjl (const RankTwoTensorTempl< T > &, const RankTwoTensorTempl< T > &, const RankTwoTensorTempl< T > &) const
	Calculates C_mnkt A_im B_jn C_lt.

RankFourTensorTempl< T >	tripleProductIjk (const RankTwoTensorTempl< T > &, const RankTwoTensorTempl< T > &, const RankTwoTensorTempl< T > &) const
	Calculates C_mntl A_im B_jn C_kt.

RankFourTensorTempl< T >	singleProductI (const RankTwoTensorTempl< T > &) const
	Calculates C_mjkl A_im.

RankFourTensorTempl< T >	singleProductJ (const RankTwoTensorTempl< T > &) const
	Calculates C_imkl A_jm.

RankFourTensorTempl< T >	singleProductK (const RankTwoTensorTempl< T > &) const
	Calculates C_ijml A_km.

RankFourTensorTempl< T >	singleProductL (const RankTwoTensorTempl< T > &) const
	Calculates C_ijkm A_lm.

bool	isSymmetric () const
	checks if the tensor is symmetric

bool	isIsotropic () const
	checks if the tensor is isotropic

RankFourTensorTempl< Real >	invSymm () const


void	fillGeneralIsotropic (const T &i0, const T &i1, const T &i2)
	Vector-less fill API functions. See docs of the corresponding ...FromInputVector methods.

void	fillAntisymmetricIsotropic (const T &i0)

void	fillSymmetricIsotropic (const T &i0, const T &i1)

void	fillSymmetricIsotropicEandNu (const T &E, const T &nu)

Static Public Member Functions
static RankFourTensorTempl< T >	Identity ()

static RankFourTensorTempl< T >	IdentityFour ()

static RankFourTensorTempl< T >	IdentityDeviatoric ()
	Identity of type \delta_{ik} \delta_{jl} - \delta_{ij} \delta_{kl} / 3.

static MooseEnum	fillMethodEnum ()
	Static method for use in validParams for getting the "fill_method".

Static Public Attributes

static constexpr unsigned int	N = Moose::dim
	tensor dimension and powers of the dimension

static constexpr unsigned int	N2 = N * N

static constexpr unsigned int	N3 = N * N * N

static constexpr unsigned int	N4 = N * N * N * N

Protected Member Functions
void	fillAntisymmetricFromInputVector (const std::vector< T > &input)
	fillAntisymmetricFromInputVector takes 6 inputs to fill the the antisymmetric Rank-4 tensor with the appropriate symmetries maintained.

void	fillGeneralIsotropicFromInputVector (const std::vector< T > &input)
	fillGeneralIsotropicFromInputVector takes 3 inputs to fill the Rank-4 tensor with symmetries C_ijkl = C_klij, and isotropy, ie C_ijkl = lade_ijde_kl + mu(de_ikde_jl + de_ilde_jk) + aep_ijm*ep_klm where la is the first Lame modulus, mu is the second (shear) Lame modulus, and a is the antisymmetric shear modulus, and ep is the permutation tensor

void	fillAntisymmetricIsotropicFromInputVector (const std::vector< T > &input)
	fillAntisymmetricIsotropicFromInputVector takes 1 input to fill the the antisymmetric Rank-4 tensor with the appropriate symmetries maintained.

void	fillSymmetricIsotropicFromInputVector (const std::vector< T > &input)
	fillSymmetricIsotropicFromInputVector takes 2 inputs to fill the the symmetric Rank-4 tensor with the appropriate symmetries maintained.

void	fillSymmetricIsotropicEandNuFromInputVector (const std::vector< T > &input)
	fillSymmetricIsotropicEandNuFromInputVector is a variation of the fillSymmetricIsotropicFromInputVector which takes as inputs the more commonly used Young's modulus (E) and Poisson's ratio (nu) constants to fill the isotropic elasticity tensor.

void	fillAxisymmetricRZFromInputVector (const std::vector< T > &input)
	fillAxisymmetricRZFromInputVector takes 5 inputs to fill the axisymmetric Rank-4 tensor with the appropriate symmetries maintatined for use with axisymmetric problems using coord_type = RZ.

void	fillGeneralFromInputVector (const std::vector< T > &input)
	fillGeneralFromInputVector takes 81 inputs to fill the Rank-4 tensor No symmetries are explicitly maintained

void	fillPrincipalFromInputVector (const std::vector< T > &input)
	fillPrincipalFromInputVector takes 9 inputs to fill a Rank-4 tensor C1111 = input0 C1122 = input1 C1133 = input2 C2211 = input3 C2222 = input4 C2233 = input5 C3311 = input6 C3322 = input7 C3333 = input8 with all other components being zero

void	fillGeneralOrthotropicFromInputVector (const std::vector< T > &input)
	fillGeneralOrhotropicFromInputVector takes 10 inputs to fill the Rank-4 tensor It defines a general orthotropic tensor for which some constraints among elastic parameters exist

Protected Attributes
T	_vals [N4]
	The values of the rank-four tensor stored by index=(((i * LIBMESH_DIM + j) * LIBMESH_DIM + k) * LIBMESH_DIM + l)

Friends
template<typename T2 >
class	RankTwoTensorTempl

template<typename T2 >
class	RankFourTensorTempl

template<typename T2 >
class	RankThreeTensorTempl

std::ostream &	operator<< (std::ostream &os, const RankFourTensorTempl< T > &t)

template<class T2 >
void	dataStore (std::ostream &, RankFourTensorTempl< T2 > &, void *)

template<class T2 >
void	dataLoad (std::istream &, RankFourTensorTempl< T2 > &, void *)

Detailed Description

template<typename T>
class RankFourTensorTempl< T >

RankFourTensorTempl is designed to handle any N-dimensional fourth order tensor, C.

Since N is hard-coded to 3, RankFourTensorTempl holds 81 separate C_ijkl entries, with i,j,k,l = 0, 1, 2.

Definition at line 65 of file RankFourTensor.h.

Member Typedef Documentation

◆ index_type

template<typename T >

typedef tuple_of<4, unsigned int> RankFourTensorTempl< T >::index_type

Definition at line 75 of file RankFourTensor.h.

◆ value_type

template<typename T >

typedef T RankFourTensorTempl< T >::value_type

Definition at line 76 of file RankFourTensor.h.

Member Enumeration Documentation

◆ FillMethod

template<typename T >

enum RankFourTensorTempl::FillMethod

To fill up the 81 entries in the 4th-order tensor, fillFromInputVector is called with one of the following fill_methods.

See the fill*FromInputVector functions for more details

Enumerator
antisymmetric
symmetric9
symmetric21
general_isotropic
symmetric_isotropic
symmetric_isotropic_E_nu
antisymmetric_isotropic
axisymmetric_rz
general
principal
orthotropic

Definition at line 93 of file RankFourTensor.h.

  {
    antisymmetric,
    symmetric9,
    symmetric21,
    general_isotropic,
    symmetric_isotropic,
    symmetric_isotropic_E_nu,
    antisymmetric_isotropic,
    axisymmetric_rz,
    general,
    principal,
    orthotropic
  };

◆ InitMethod

template<typename T >

enum RankFourTensorTempl::InitMethod

Initialization method.

Enumerator
initNone
initIdentity
initIdentityFour
initIdentitySymmetricFour
initIdentityDeviatoric

Definition at line 79 of file RankFourTensor.h.

  {
    initNone,
    initIdentity,
    initIdentityFour,
    initIdentitySymmetricFour,
    initIdentityDeviatoric
  };

Constructor & Destructor Documentation

◆ RankFourTensorTempl() [1/6]

template<typename T >

RankFourTensorTempl< T >::RankFourTensorTempl ( )

Default constructor; fills to zero.

Definition at line 52 of file RankFourTensorImplementation.h.

{
  mooseAssert(N == 3, "RankFourTensorTempl<T> is currently only tested for 3 dimensions.");
 
  for (auto i : make_range(N4))
    _vals[i] = 0.0;
}

◆ RankFourTensorTempl() [2/6]

template<typename T >

RankFourTensorTempl< T >::RankFourTensorTempl ( const InitMethod init )

Select specific initialization pattern.

Definition at line 61 of file RankFourTensorImplementation.h.

{
  unsigned int index = 0;
  switch (init)
  {
    case initNone:
      break;
 
    case initIdentity:
      zero();
      for (auto i : make_range(N))
        (*this)(i, i, i, i) = 1.0;
      break;
 
    case initIdentityFour:
      for (auto i : make_range(N))
        for (auto j : make_range(N))
          for (auto k : make_range(N))
            for (auto l : make_range(N))
              _vals[index++] = Real(i == k && j == l);
      break;
 
    case initIdentitySymmetricFour:
      for (auto i : make_range(N))
        for (auto j : make_range(N))
          for (auto k : make_range(N))
            for (auto l : make_range(N))
              _vals[index++] = 0.5 * Real(i == k && j == l) + 0.5 * Real(i == l && j == k);
      break;
 
    case initIdentityDeviatoric:
      for (unsigned int i = 0; i < N; ++i)
        for (unsigned int j = 0; j < N; ++j)
          for (unsigned int k = 0; k < N; ++k)
            for (unsigned int l = 0; l < N; ++l)
            {
              _vals[index] = Real(i == k && j == l);
              if ((i == j) && (k == l))
                _vals[index] -= 1.0 / 3.0;
              index++;
            }
      break;
 
    default:
      mooseError("Unknown RankFourTensorTempl<T> initialization pattern.");
  }
}

◆ RankFourTensorTempl() [3/6]

template<typename T >

RankFourTensorTempl< T >::RankFourTensorTempl	(	const std::vector< T > &	input,
		FillMethod	fill_method
	)

Fill from vector.

Definition at line 110 of file RankFourTensorImplementation.h.

{
  fillFromInputVector(input, fill_method);
}

◆ RankFourTensorTempl() [4/6]

template<typename T >

RankFourTensorTempl< T >::RankFourTensorTempl ( const RankFourTensorTempl< T > & a )

default

Copy assignment operator must be defined if used.

◆ RankFourTensorTempl() [5/6]

template<typename T >

template<typename T2 >

RankFourTensorTempl< T >::RankFourTensorTempl ( const RankFourTensorTempl< T2 > & copy )

Copy constructor.

Definition at line 570 of file RankFourTensor.h.

{
  for (auto i : libMesh::make_range(N4))
    _vals[i] = copy._vals[i];
}

◆ RankFourTensorTempl() [6/6]

template<typename T >

template<typename T2 >

RankFourTensorTempl< T >::RankFourTensorTempl ( const SymmetricRankFourTensorTempl< T2 > & t )

inline

The conversion operator from a SymmetricRankFourTensorTempl

Definition at line 151 of file RankFourTensor.h.

  {
    for (const auto a : make_range(SymmetricRankFourTensorTempl<T2>::N))
      for (const auto b : make_range(SymmetricRankFourTensorTempl<T2>::N))
      {
        const auto & idx = SymmetricRankFourTensorTempl<T2>::full_index[a][b];
        const auto i = idx[0];
        const auto j = idx[1];
        const auto k = idx[2];
        const auto l = idx[3];
 
        // Rijkl = Rjikl = Rijlk = Rjilk
        (*this)(i, j, k, l) = t(a, b) / SymmetricRankFourTensorTempl<T2>::mandelFactor(a, b);
        (*this)(j, i, k, l) = t(a, b) / SymmetricRankFourTensorTempl<T2>::mandelFactor(a, b);
        (*this)(i, j, l, k) = t(a, b) / SymmetricRankFourTensorTempl<T2>::mandelFactor(a, b);
        (*this)(j, i, l, k) = t(a, b) / SymmetricRankFourTensorTempl<T2>::mandelFactor(a, b);
      }
  }

Member Function Documentation

◆ contraction()

template<typename T >

template<int m>

RankThreeTensorTempl< T > RankFourTensorTempl< T >::contraction ( const libMesh::VectorValue< T > & b ) const

single contraction of a RankFourTensor with a vector over index m

Returns: C_xxx = a_ijkl*b_m where m={i,j,k,l} and xxx the remaining indices

Definition at line 789 of file RankFourTensor.h.

{
  RankThreeTensorTempl<T> result;
  static constexpr std::size_t z[4][3] = {{1, 2, 3}, {0, 2, 3}, {0, 1, 3}, {0, 1, 2}};
  std::size_t x[4];
  for (x[0] = 0; x[0] < N; ++x[0])
    for (x[1] = 0; x[1] < N; ++x[1])
      for (x[2] = 0; x[2] < N; ++x[2])
        for (x[3] = 0; x[3] < N; ++x[3])
          result(x[z[m][0]], x[z[m][1]], x[z[m][2]]) += (*this)(x[0], x[1], x[2], x[3]) * b(x[m]);
 
  return result;
}

◆ contractionIj()

template<typename T >

T RankFourTensorTempl< T >::contractionIj	(	unsigned int	i,
		unsigned int	j,
		const RankTwoTensorTempl< T > &	M
	)		const

Sum C_ijkl M_kl for a given i,j.

Definition at line 917 of file RankFourTensorImplementation.h.

{
  T val = 0;
  for (unsigned int k = 0; k < N; k++)
    for (unsigned int l = 0; l < N; l++)
      val += (*this)(i, j, k, l) * M(k, l);
 
  return val;
}

◆ contractionKl()

template<typename T >

T RankFourTensorTempl< T >::contractionKl	(	unsigned int	k,
		unsigned int	l,
		const RankTwoTensorTempl< T > &	M
	)		const

Sum M_ij C_ijkl for a given k,l.

Definition at line 931 of file RankFourTensorImplementation.h.

{
  T val = 0;
  for (unsigned int i = 0; i < N; i++)
    for (unsigned int j = 0; j < N; j++)
      val += (*this)(i, j, k, l) * M(i, j);
 
  return val;
}

◆ fillAntisymmetricFromInputVector()

template<typename T >

void RankFourTensorTempl< T >::fillAntisymmetricFromInputVector ( const std::vector< T > & input )

protected

fillAntisymmetricFromInputVector takes 6 inputs to fill the the antisymmetric Rank-4 tensor with the appropriate symmetries maintained.

I.e., B_ijkl = -B_jikl = -B_ijlk = B_klij

Parameters

input this is B1212, B1213, B1223, B1313, B1323, B2323.

Definition at line 621 of file RankFourTensorImplementation.h.

{
  if (input.size() != 6)
    mooseError(
        "To use fillAntisymmetricFromInputVector, your input must have size 6.  Yours has size ",
        input.size());
 
  zero();
 
  (*this)(0, 1, 0, 1) = input[0]; // B1212
  (*this)(0, 1, 0, 2) = input[1]; // B1213
  (*this)(0, 1, 1, 2) = input[2]; // B1223
 
  (*this)(0, 2, 0, 2) = input[3]; // B1313
  (*this)(0, 2, 1, 2) = input[4]; // B1323
 
  (*this)(1, 2, 1, 2) = input[5]; // B2323
 
  // symmetry on the two pairs
  (*this)(0, 2, 0, 1) = (*this)(0, 1, 0, 2);
  (*this)(1, 2, 0, 1) = (*this)(0, 1, 1, 2);
  (*this)(1, 2, 0, 2) = (*this)(0, 2, 1, 2);
  // have now got the upper parts of vals[0][1], vals[0][2] and vals[1][2]
 
  // fill in from antisymmetry relations
  for (auto i : make_range(N))
    for (auto j : make_range(N))
    {
      (*this)(0, 1, j, i) = -(*this)(0, 1, i, j);
      (*this)(0, 2, j, i) = -(*this)(0, 2, i, j);
      (*this)(1, 2, j, i) = -(*this)(1, 2, i, j);
    }
  // have now got all of vals[0][1], vals[0][2] and vals[1][2]
 
  // fill in from antisymmetry relations
  for (auto i : make_range(N))
    for (auto j : make_range(N))
    {
      (*this)(1, 0, i, j) = -(*this)(0, 1, i, j);
      (*this)(2, 0, i, j) = -(*this)(0, 2, i, j);
      (*this)(2, 1, i, j) = -(*this)(1, 2, i, j);
    }
}

◆ fillAntisymmetricIsotropic()

template<typename T >

void RankFourTensorTempl< T >::fillAntisymmetricIsotropic ( const T & i0 )

Definition at line 708 of file RankFourTensorImplementation.h.

{
  fillGeneralIsotropic(0.0, 0.0, i0);
}

◆ fillAntisymmetricIsotropicFromInputVector()

template<typename T >

void RankFourTensorTempl< T >::fillAntisymmetricIsotropicFromInputVector ( const std::vector< T > & input )

protected

fillAntisymmetricIsotropicFromInputVector takes 1 input to fill the the antisymmetric Rank-4 tensor with the appropriate symmetries maintained.

I.e., C_ijkl = a * ep_ijm * ep_klm, where epsilon is the permutation tensor (and sum on m)

Parameters

input this is a in the above formula

Definition at line 696 of file RankFourTensorImplementation.h.

{
  if (input.size() != 1)
    mooseError("To use fillAntisymmetricIsotropicFromInputVector, your input must have size 1. "
               "Yours has size ",
               input.size());
 
  fillGeneralIsotropic(0.0, 0.0, input[0]);
}

◆ fillAxisymmetricRZFromInputVector()

template<typename T >

void RankFourTensorTempl< T >::fillAxisymmetricRZFromInputVector ( const std::vector< T > & input )

protected

fillAxisymmetricRZFromInputVector takes 5 inputs to fill the axisymmetric Rank-4 tensor with the appropriate symmetries maintatined for use with axisymmetric problems using coord_type = RZ.

I.e. C1111 = C2222, C1133 = C2233, C2323 = C3131 and C1212 = 0.5*(C1111-C1122)

Parameters

input this is C1111, C1122, C1133, C3333, C2323.

Definition at line 763 of file RankFourTensorImplementation.h.

{
  mooseAssert(input.size() == 5,
              "To use fillAxisymmetricRZFromInputVector, your input must have size 5.");
 
  // C1111  C1122     C1133       0         0         0
  //        C2222  C2233=C1133    0         0         0
  //                  C3333       0         0         0
  //                            C2323       0         0
  //                                   C3131=C2323    0
  //                                                C1212
  // clang-format off
  fillSymmetric21FromInputVector(std::array<T,21>
  {{input[0],input[1],input[2],      0.0,      0.0, 0.0,
             input[0],input[2],      0.0,      0.0, 0.0,
                      input[3],      0.0,      0.0, 0.0,
                                input[4],      0.0, 0.0,
                                          input[4], 0.0,
                                                    (input[0] - input[1]) * 0.5}});
  // clang-format on
}

◆ fillFromInputVector()

template<typename T >

void RankFourTensorTempl< T >::fillFromInputVector	(	const std::vector< T > &	input,
		FillMethod	fill_method
	)

fillFromInputVector takes some number of inputs to fill the Rank-4 tensor.

Parameters

input the numbers that will be placed in the tensor

fill_method this can be: antisymmetric (use fillAntisymmetricFromInputVector) symmetric9 (use fillSymmetric9FromInputVector) symmetric21 (use fillSymmetric21FromInputVector) general_isotropic (use fillGeneralIsotropicFrominputVector) symmetric_isotropic (use fillSymmetricIsotropicFromInputVector) antisymmetric_isotropic (use fillAntisymmetricIsotropicFromInputVector) axisymmetric_rz (use fillAxisymmetricRZFromInputVector) general (use fillGeneralFromInputVector) principal (use fillPrincipalFromInputVector)

Definition at line 576 of file RankFourTensorImplementation.h.

{
 
  switch (fill_method)
  {
    case antisymmetric:
      fillAntisymmetricFromInputVector(input);
      break;
    case symmetric9:
      fillSymmetric9FromInputVector(input);
      break;
    case symmetric21:
      fillSymmetric21FromInputVector(input);
      break;
    case general_isotropic:
      fillGeneralIsotropicFromInputVector(input);
      break;
    case symmetric_isotropic:
      fillSymmetricIsotropicFromInputVector(input);
      break;
    case symmetric_isotropic_E_nu:
      fillSymmetricIsotropicEandNuFromInputVector(input);
      break;
    case antisymmetric_isotropic:
      fillAntisymmetricIsotropicFromInputVector(input);
      break;
    case axisymmetric_rz:
      fillAxisymmetricRZFromInputVector(input);
      break;
    case general:
      fillGeneralFromInputVector(input);
      break;
    case principal:
      fillPrincipalFromInputVector(input);
      break;
    case orthotropic:
      fillGeneralOrthotropicFromInputVector(input);
      break;
    default:
      mooseError("fillFromInputVector called with unknown fill_method of ", fill_method);
  }
}

◆ fillGeneralFromInputVector()

template<typename T >

void RankFourTensorTempl< T >::fillGeneralFromInputVector ( const std::vector< T > & input )

protected

fillGeneralFromInputVector takes 81 inputs to fill the Rank-4 tensor No symmetries are explicitly maintained

Parameters

input C(i,j,k,l) = input[i*N*N*N + j*N*N + k*N + l]

Definition at line 787 of file RankFourTensorImplementation.h.

{
  if (input.size() != 81)
    mooseError("To use fillGeneralFromInputVector, your input must have size 81. Yours has size ",
               input.size());
 
  for (auto i : make_range(N4))
    _vals[i] = input[i];
}

◆ fillGeneralIsotropic()

template<typename T >

void RankFourTensorTempl< T >::fillGeneralIsotropic	(	const T &	i0,
		const T &	i1,
		const T &	i2
	)

Vector-less fill API functions. See docs of the corresponding ...FromInputVector methods.

Definition at line 679 of file RankFourTensorImplementation.h.

{
  for (auto i : make_range(N))
    for (auto j : make_range(N))
      for (auto k : make_range(N))
        for (auto l : make_range(N))
        {
          (*this)(i, j, k, l) = i0 * Real(i == j) * Real(k == l) +
                                i1 * Real(i == k) * Real(j == l) + i1 * Real(i == l) * Real(j == k);
          for (auto m : make_range(N))
            (*this)(i, j, k, l) +=
                i2 * Real(PermutationTensor::eps(i, j, m)) * Real(PermutationTensor::eps(k, l, m));
        }
}

◆ fillGeneralIsotropicFromInputVector()

template<typename T >

void RankFourTensorTempl< T >::fillGeneralIsotropicFromInputVector ( const std::vector< T > & input )

protected

fillGeneralIsotropicFromInputVector takes 3 inputs to fill the Rank-4 tensor with symmetries C_ijkl = C_klij, and isotropy, ie C_ijkl = la*de_ij*de_kl + mu*(de_ik*de_jl + de_il*de_jk) + a*ep_ijm*ep_klm where la is the first Lame modulus, mu is the second (shear) Lame modulus, and a is the antisymmetric shear modulus, and ep is the permutation tensor

Parameters

input this is la, mu, a in the above formula

Definition at line 667 of file RankFourTensorImplementation.h.

{
  if (input.size() != 3)
    mooseError("To use fillGeneralIsotropicFromInputVector, your input must have size 3.  Yours "
               "has size ",
               input.size());
 
  fillGeneralIsotropic(input[0], input[1], input[2]);
}

◆ fillGeneralOrthotropicFromInputVector()

template<typename T >

void RankFourTensorTempl< T >::fillGeneralOrthotropicFromInputVector ( const std::vector< T > & input )

protected

fillGeneralOrhotropicFromInputVector takes 10 inputs to fill the Rank-4 tensor It defines a general orthotropic tensor for which some constraints among elastic parameters exist

Parameters

input Ea, Eb, Ec, Gab, Gbc, Gca, nuba, nuca, nucb, nuab, nuac, nubc

Definition at line 820 of file RankFourTensorImplementation.h.

{
  if (input.size() != 12)
    mooseError("To use fillGeneralOrhotropicFromInputVector, your input must have size 12. Yours "
               "has size ",
               input.size());
 
  const T & Ea = input[0];
  const T & Eb = input[1];
  const T & Ec = input[2];
  const T & Gab = input[3];
  const T & Gbc = input[4];
  const T & Gca = input[5];
  const T & nuba = input[6];
  const T & nuca = input[7];
  const T & nucb = input[8];
  const T & nuab = input[9];
  const T & nuac = input[10];
  const T & nubc = input[11];
 
  // Input must satisfy constraints.
  bool preserve_symmetry = MooseUtils::relativeFuzzyEqual(nuab * Eb, nuba * Ea) &&
                           MooseUtils::relativeFuzzyEqual(nuca * Ea, nuac * Ec) &&
                           MooseUtils::relativeFuzzyEqual(nubc * Ec, nucb * Eb);
 
  if (!preserve_symmetry)
    mooseError("Orthotropic elasticity tensor input is not consistent with symmetry requirements. "
               "Check input for accuracy");
 
  unsigned int ntens = N * (N + 1) / 2;
 
  std::vector<T> mat;
  mat.assign(ntens * ntens, 0);
 
  T k = 1 - nuab * nuba - nubc * nucb - nuca * nuac - nuab * nubc * nuca - nuba * nucb * nuac;
 
  bool is_positive_definite =
      (k > 0) && (1 - nubc * nucb) > 0 && (1 - nuac * nuca) > 0 && (1 - nuab * nuba) > 0;
  if (!is_positive_definite)
    mooseError("Orthotropic elasticity tensor input is not positive definite. Check input for "
               "accuracy");
 
  mat[0] = Ea * (1 - nubc * nucb) / k;
  mat[1] = Ea * (nubc * nuca + nuba) / k;
  mat[2] = Ea * (nuba * nucb + nuca) / k;
 
  mat[6] = Eb * (nuac * nucb + nuab) / k;
  mat[7] = Eb * (1 - nuac * nuca) / k;
  mat[8] = Eb * (nuab * nuca + nucb) / k;
 
  mat[12] = Ec * (nuab * nubc + nuac) / k;
  mat[13] = Ec * (nuac * nuba + nubc) / k;
  mat[14] = Ec * (1 - nuab * nuba) / k;
 
  mat[21] = 2 * Gab;
  mat[28] = 2 * Gca;
  mat[35] = 2 * Gbc;
 
  // Switching from Voigt to fourth order tensor
  // Copied from existing code (invSymm)
  int nskip = N - 1;
 
  unsigned int index = 0;
  for (auto i : make_range(N))
    for (auto j : make_range(N))
      for (auto k : make_range(N))
        for (auto l : make_range(N))
        {
          if (i == j)
            (*this)._vals[index] =
                k == l ? mat[i * ntens + k] : mat[i * ntens + k + nskip + l] / 2.0;
          else
            (*this)._vals[index] = k == l ? mat[(nskip + i + j) * ntens + k]
                                          : mat[(nskip + i + j) * ntens + k + nskip + l] / 2.0;
          index++;
        }
}

◆ fillMethodEnum()

template<typename T >

MooseEnum RankFourTensorTempl< T >::fillMethodEnum ( )

static

Static method for use in validParams for getting the "fill_method".

Definition at line 44 of file RankFourTensorImplementation.h.

{
  return MooseEnum("antisymmetric symmetric9 symmetric21 general_isotropic symmetric_isotropic "
                   "symmetric_isotropic_E_nu antisymmetric_isotropic axisymmetric_rz general "
                   "principal orthotropic");
}

◆ fillPrincipalFromInputVector()

template<typename T >

void RankFourTensorTempl< T >::fillPrincipalFromInputVector ( const std::vector< T > & input )

protected

fillPrincipalFromInputVector takes 9 inputs to fill a Rank-4 tensor C1111 = input0 C1122 = input1 C1133 = input2 C2211 = input3 C2222 = input4 C2233 = input5 C3311 = input6 C3322 = input7 C3333 = input8 with all other components being zero

Definition at line 799 of file RankFourTensorImplementation.h.

{
  if (input.size() != 9)
    mooseError("To use fillPrincipalFromInputVector, your input must have size 9. Yours has size ",
               input.size());
 
  zero();
 
  (*this)(0, 0, 0, 0) = input[0];
  (*this)(0, 0, 1, 1) = input[1];
  (*this)(0, 0, 2, 2) = input[2];
  (*this)(1, 1, 0, 0) = input[3];
  (*this)(1, 1, 1, 1) = input[4];
  (*this)(1, 1, 2, 2) = input[5];
  (*this)(2, 2, 0, 0) = input[6];
  (*this)(2, 2, 1, 1) = input[7];
  (*this)(2, 2, 2, 2) = input[8];
}

◆ fillSymmetric21FromInputVector()

template<typename T >

template<typename T2 >

void RankFourTensorTempl< T >::fillSymmetric21FromInputVector ( const T2 & input )

fillSymmetric21FromInputVector takes either 21 inputs to fill in the Rank-4 tensor with the appropriate crystal symmetries maintained.

I.e., C_ijkl = C_klij, C_ijkl = C_ijlk, C_ijkl = C_jikl

Parameters

input is C1111 C1122 C1133 C1123 C1113 C1112 C2222 C2233 C2223 C2213 C2212 C3333 C3323 C3313 C3312 C2323 C2313 C2312 C1313 C1312 C1212

Definition at line 645 of file RankFourTensor.h.

{
  mooseAssert(input.size() == 21,
              "To use fillSymmetric21FromInputVector, your input must have size 21.");
 
  (*this)(0, 0, 0, 0) = input[0];  // C1111
  (*this)(1, 1, 1, 1) = input[6];  // C2222
  (*this)(2, 2, 2, 2) = input[11]; // C3333
 
  (*this)(0, 0, 1, 1) = input[1]; // C1122
  (*this)(1, 1, 0, 0) = input[1];
 
  (*this)(0, 0, 2, 2) = input[2]; // C1133
  (*this)(2, 2, 0, 0) = input[2];
 
  (*this)(1, 1, 2, 2) = input[7]; // C2233
  (*this)(2, 2, 1, 1) = input[7];
 
  (*this)(0, 0, 0, 2) = input[4]; // C1113
  (*this)(0, 0, 2, 0) = input[4];
  (*this)(0, 2, 0, 0) = input[4];
  (*this)(2, 0, 0, 0) = input[4];
 
  (*this)(0, 0, 0, 1) = input[5]; // C1112
  (*this)(0, 0, 1, 0) = input[5];
  (*this)(0, 1, 0, 0) = input[5];
  (*this)(1, 0, 0, 0) = input[5];
 
  (*this)(1, 1, 1, 2) = input[8]; // C2223
  (*this)(1, 1, 2, 1) = input[8];
  (*this)(1, 2, 1, 1) = input[8];
  (*this)(2, 1, 1, 1) = input[8];
 
  (*this)(1, 1, 1, 0) = input[10];
  (*this)(1, 1, 0, 1) = input[10];
  (*this)(1, 0, 1, 1) = input[10];
  (*this)(0, 1, 1, 1) = input[10]; // C2212 //flipped for filling purposes
 
  (*this)(2, 2, 2, 1) = input[12];
  (*this)(2, 2, 1, 2) = input[12];
  (*this)(2, 1, 2, 2) = input[12];
  (*this)(1, 2, 2, 2) = input[12]; // C3323 //flipped for filling purposes
 
  (*this)(2, 2, 2, 0) = input[13];
  (*this)(2, 2, 0, 2) = input[13];
  (*this)(2, 0, 2, 2) = input[13];
  (*this)(0, 2, 2, 2) = input[13]; // C3313 //flipped for filling purposes
 
  (*this)(0, 0, 1, 2) = input[3]; // C1123
  (*this)(0, 0, 2, 1) = input[3];
  (*this)(1, 2, 0, 0) = input[3];
  (*this)(2, 1, 0, 0) = input[3];
 
  (*this)(1, 1, 0, 2) = input[9];
  (*this)(1, 1, 2, 0) = input[9];
  (*this)(0, 2, 1, 1) = input[9]; // C2213  //flipped for filling purposes
  (*this)(2, 0, 1, 1) = input[9];
 
  (*this)(2, 2, 0, 1) = input[14];
  (*this)(2, 2, 1, 0) = input[14];
  (*this)(0, 1, 2, 2) = input[14]; // C3312 //flipped for filling purposes
  (*this)(1, 0, 2, 2) = input[14];
 
  (*this)(1, 2, 1, 2) = input[15]; // C2323
  (*this)(2, 1, 2, 1) = input[15];
  (*this)(2, 1, 1, 2) = input[15];
  (*this)(1, 2, 2, 1) = input[15];
 
  (*this)(0, 2, 0, 2) = input[18]; // C1313
  (*this)(2, 0, 2, 0) = input[18];
  (*this)(2, 0, 0, 2) = input[18];
  (*this)(0, 2, 2, 0) = input[18];
 
  (*this)(0, 1, 0, 1) = input[20]; // C1212
  (*this)(1, 0, 1, 0) = input[20];
  (*this)(1, 0, 0, 1) = input[20];
  (*this)(0, 1, 1, 0) = input[20];
 
  (*this)(1, 2, 0, 2) = input[16];
  (*this)(0, 2, 1, 2) = input[16]; // C2313 //flipped for filling purposes
  (*this)(2, 1, 0, 2) = input[16];
  (*this)(1, 2, 2, 0) = input[16];
  (*this)(2, 0, 1, 2) = input[16];
  (*this)(0, 2, 2, 1) = input[16];
  (*this)(2, 1, 2, 0) = input[16];
  (*this)(2, 0, 2, 1) = input[16];
 
  (*this)(1, 2, 0, 1) = input[17];
  (*this)(0, 1, 1, 2) = input[17]; // C2312 //flipped for filling purposes
  (*this)(2, 1, 0, 1) = input[17];
  (*this)(1, 2, 1, 0) = input[17];
  (*this)(1, 0, 1, 2) = input[17];
  (*this)(0, 1, 2, 1) = input[17];
  (*this)(2, 1, 1, 0) = input[17];
  (*this)(1, 0, 2, 1) = input[17];
 
  (*this)(0, 2, 0, 1) = input[19];
  (*this)(0, 1, 0, 2) = input[19]; // C1312 //flipped for filling purposes
  (*this)(2, 0, 0, 1) = input[19];
  (*this)(0, 2, 1, 0) = input[19];
  (*this)(1, 0, 0, 2) = input[19];
  (*this)(0, 1, 2, 0) = input[19];
  (*this)(2, 0, 1, 0) = input[19];
  (*this)(1, 0, 2, 0) = input[19];
}

◆ fillSymmetric9FromInputVector()

template<typename T >

template<typename T2 >

void RankFourTensorTempl< T >::fillSymmetric9FromInputVector ( const T2 & input )

fillSymmetric9FromInputVector takes 9 inputs to fill in the Rank-4 tensor with the appropriate crystal symmetries maintained.

I.e., C_ijkl = C_klij, C_ijkl = C_ijlk, C_ijkl = C_jikl

Parameters

input is: C1111 C1122 C1133 C2222 C2233 C3333 C2323 C1313 C1212 In the isotropic case this is (la is first Lame constant, mu is second (shear) Lame constant) la+2mu la la la+2mu la la+2mu mu mu mu

Definition at line 608 of file RankFourTensor.h.

{
  mooseAssert(input.size() == 9,
              "To use fillSymmetric9FromInputVector, your input must have size 9.");
  zero();
 
  (*this)(0, 0, 0, 0) = input[0]; // C1111
  (*this)(1, 1, 1, 1) = input[3]; // C2222
  (*this)(2, 2, 2, 2) = input[5]; // C3333
 
  (*this)(0, 0, 1, 1) = input[1]; // C1122
  (*this)(1, 1, 0, 0) = input[1];
 
  (*this)(0, 0, 2, 2) = input[2]; // C1133
  (*this)(2, 2, 0, 0) = input[2];
 
  (*this)(1, 1, 2, 2) = input[4]; // C2233
  (*this)(2, 2, 1, 1) = input[4];
 
  (*this)(1, 2, 1, 2) = input[6]; // C2323
  (*this)(2, 1, 2, 1) = input[6];
  (*this)(2, 1, 1, 2) = input[6];
  (*this)(1, 2, 2, 1) = input[6];
 
  (*this)(0, 2, 0, 2) = input[7]; // C1313
  (*this)(2, 0, 2, 0) = input[7];
  (*this)(2, 0, 0, 2) = input[7];
  (*this)(0, 2, 2, 0) = input[7];
 
  (*this)(0, 1, 0, 1) = input[8]; // C1212
  (*this)(1, 0, 1, 0) = input[8];
  (*this)(1, 0, 0, 1) = input[8];
  (*this)(0, 1, 1, 0) = input[8];
}

◆ fillSymmetricIsotropic()

template<typename T >

void RankFourTensorTempl< T >::fillSymmetricIsotropic	(	const T &	i0,
		const T &	i1
	)

Definition at line 724 of file RankFourTensorImplementation.h.

{
  // clang-format off
  fillSymmetric21FromInputVector(std::array<T,21>
  {{lambda + 2.0 * G, lambda,           lambda,           0.0, 0.0, 0.0,
                      lambda + 2.0 * G, lambda,           0.0, 0.0, 0.0,
                                        lambda + 2.0 * G, 0.0, 0.0, 0.0,
                                                            G, 0.0, 0.0,
                                                                 G, 0.0,
                                                                      G}});
  // clang-format on
}

◆ fillSymmetricIsotropicEandNu()

template<typename T >

void RankFourTensorTempl< T >::fillSymmetricIsotropicEandNu	(	const T &	E,
		const T &	nu
	)

Definition at line 752 of file RankFourTensorImplementation.h.

{
  // Calculate lambda and the shear modulus from the given young's modulus and poisson's ratio
  const T & lambda = E * nu / ((1.0 + nu) * (1.0 - 2.0 * nu));
  const T & G = E / (2.0 * (1.0 + nu));
 
  fillSymmetricIsotropic(lambda, G);
}

◆ fillSymmetricIsotropicEandNuFromInputVector()

template<typename T >

void RankFourTensorTempl< T >::fillSymmetricIsotropicEandNuFromInputVector ( const std::vector< T > & input )

protected

fillSymmetricIsotropicEandNuFromInputVector is a variation of the fillSymmetricIsotropicFromInputVector which takes as inputs the more commonly used Young's modulus (E) and Poisson's ratio (nu) constants to fill the isotropic elasticity tensor.

Using well-known formulas, E and nu are used to calculate lambda and mu and then the vector is passed to fillSymmetricIsotropicFromInputVector.

Parameters

input Young's modulus (E) and Poisson's ratio (nu)

Definition at line 739 of file RankFourTensorImplementation.h.

{
  if (input.size() != 2)
    mooseError(
        "To use fillSymmetricIsotropicEandNuFromInputVector, your input must have size 2. Yours "
        "has size ",
        input.size());
 
  fillSymmetricIsotropicEandNu(input[0], input[1]);
}

◆ fillSymmetricIsotropicFromInputVector()

template<typename T >

void RankFourTensorTempl< T >::fillSymmetricIsotropicFromInputVector ( const std::vector< T > & input )

protected

fillSymmetricIsotropicFromInputVector takes 2 inputs to fill the the symmetric Rank-4 tensor with the appropriate symmetries maintained.

C_ijkl = lambda*de_ij*de_kl + mu*(de_ik*de_jl + de_il*de_jk) where lambda is the first Lame modulus, mu is the second (shear) Lame modulus,

Parameters

input this is lambda and mu in the above formula

Definition at line 715 of file RankFourTensorImplementation.h.

{
  mooseAssert(input.size() == 2,
              "To use fillSymmetricIsotropicFromInputVector, your input must have size 2.");
  fillSymmetricIsotropic(input[0], input[1]);
}

◆ Identity()

template<typename T >

static RankFourTensorTempl< T > RankFourTensorTempl< T >::Identity ( )

inlinestatic

Definition at line 171 of file RankFourTensor.h.

171{ return RankFourTensorTempl<T>(initIdentity); }

RankFourTensorTempl

RankFourTensorTempl is designed to handle any N-dimensional fourth order tensor, C.

Definition RankFourTensor.h:66

◆ IdentityDeviatoric()

template<typename T >

static RankFourTensorTempl< T > RankFourTensorTempl< T >::IdentityDeviatoric ( )

inlinestatic

Identity of type \delta_{ik} \delta_{jl} - \delta_{ij} \delta_{kl} / 3.

Definition at line 174 of file RankFourTensor.h.

  {
    return RankFourTensorTempl<T>(initIdentityDeviatoric);
  };

◆ IdentityFour()

template<typename T >

static RankFourTensorTempl< T > RankFourTensorTempl< T >::IdentityFour ( )

inlinestatic

Definition at line 172 of file RankFourTensor.h.

172{ return RankFourTensorTempl<T>(initIdentityFour); };

◆ innerProductTranspose()

template<typename T >

RankTwoTensorTempl< T > RankFourTensorTempl< T >::innerProductTranspose ( const RankTwoTensorTempl< T > & b ) const

Inner product of the major transposed tensor with a rank two tensor.

Definition at line 900 of file RankFourTensorImplementation.h.

{
  RankTwoTensorTempl<T> result;
 
  unsigned int index = 0;
  for (unsigned int ij = 0; ij < N2; ++ij)
  {
    T bb = b._coords[ij];
    for (unsigned int kl = 0; kl < N2; ++kl)
      result._coords[kl] += _vals[index++] * bb;
  }
 
  return result;
}

◆ inverse()

template<typename T >

RankFourTensorTempl< T > RankFourTensorTempl< T >::inverse ( ) const

This returns A_ijkl such that C_ijkl*A_klmn = de_im de_jn i.e.

the general rank four inverse

Definition at line 753 of file RankFourTensor.h.

{
 
  // The inverse of a 3x3x3x3 in the C_ijkl*A_klmn = de_im de_jn sense is
  // simply the inverse of the 9x9 matrix of the tensor entries.
  // So all we need to do is inverse _vals (with the appropriate row-major
  // storage)
 
  RankFourTensorTempl<T> result;
 
  if constexpr (RankFourTensorTempl<T>::N4 * sizeof(T) > EIGEN_STACK_ALLOCATION_LIMIT)
  {
    // Allocate on the heap if you're going to exceed the stack size limit
    Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> mat(9, 9);
    for (auto i : libMesh::make_range(9 * 9))
      mat(i) = _vals[i];
 
    mat = mat.inverse();
 
    for (auto i : libMesh::make_range(9 * 9))
      result._vals[i] = mat(i);
  }
  else
  {
    // Allocate on the stack if small enough
    const Eigen::Map<const Eigen::Matrix<T, 9, 9, Eigen::RowMajor>> mat(&_vals[0]);
    Eigen::Map<Eigen::Matrix<T, 9, 9, Eigen::RowMajor>> res(&result._vals[0]);
    res = mat.inverse();
  }
 
  return result;
}

◆ invSymm() [1/2]

template<typename T >

RankFourTensorTempl< T > RankFourTensorTempl< T >::invSymm ( ) const

This returns A_ijkl such that C_ijkl*A_klmn = 0.5*(de_im de_jn + de_in de_jm) This routine assumes that C_ijkl = C_jikl = C_ijlk.

Definition at line 281 of file RankFourTensorImplementation.h.

{
  mooseError("The invSymm operation calls to LAPACK and only supports plain Real type tensors.");
}

◆ invSymm() [2/2]

RankFourTensorTempl< Real > RankFourTensorTempl< Real >::invSymm ( ) const

Definition at line 288 of file RankFourTensorImplementation.h.

{
  unsigned int ntens = N * (N + 1) / 2;
  int nskip = N - 1;
 
  RankFourTensorTempl<Real> result;
  std::vector<PetscScalar> mat;
  mat.assign(ntens * ntens, 0);
 
  // We use the LAPACK matrix inversion routine here.  Form the matrix
  //
  // mat[0]  mat[1]  mat[2]  mat[3]  mat[4]  mat[5]
  // mat[6]  mat[7]  mat[8]  mat[9]  mat[10] mat[11]
  // mat[12] mat[13] mat[14] mat[15] mat[16] mat[17]
  // mat[18] mat[19] mat[20] mat[21] mat[22] mat[23]
  // mat[24] mat[25] mat[26] mat[27] mat[28] mat[29]
  // mat[30] mat[31] mat[32] mat[33] mat[34] mat[35]
  //
  // This is filled from the indpendent components of C assuming
  // the symmetry C_ijkl = C_ijlk = C_jikl.
  //
  // If there are two rank-four tensors X and Y then the reason for
  // this filling becomes apparent if we want to calculate
  // X_ijkl*Y_klmn = Z_ijmn
  // For denote the "mat" versions of X, Y and Z by x, y and z.
  // Then
  // z_ab = x_ac*y_cb
  // Eg
  // z_00 = Z_0000 = X_0000*Y_0000 + X_0011*Y_1111 + X_0022*Y_2200 + 2*X_0001*Y_0100 +
  // 2*X_0002*Y_0200 + 2*X_0012*Y_1200   (the factors of 2 come from the assumed symmetries)
  // z_03 = 2*Z_0001 = X_0000*2*Y_0001 + X_0011*2*Y_1101 + X_0022*2*Y_2201 + 2*X_0001*2*Y_0101 +
  // 2*X_0002*2*Y_0201 + 2*X_0012*2*Y_1201
  // z_22 = 2*Z_0102 = X_0100*2*Y_0002 + X_0111*2*X_1102 + X_0122*2*Y_2202 + 2*X_0101*2*Y_0102 +
  // 2*X_0102*2*Y_0202 + 2*X_0112*2*Y_1202
  // Finally, we use LAPACK to find x^-1, and put it back into rank-4 tensor form
  //
  // mat[0] = C(0,0,0,0)
  // mat[1] = C(0,0,1,1)
  // mat[2] = C(0,0,2,2)
  // mat[3] = C(0,0,0,1)*2
  // mat[4] = C(0,0,0,2)*2
  // mat[5] = C(0,0,1,2)*2
 
  // mat[6] = C(1,1,0,0)
  // mat[7] = C(1,1,1,1)
  // mat[8] = C(1,1,2,2)
  // mat[9] = C(1,1,0,1)*2
  // mat[10] = C(1,1,0,2)*2
  // mat[11] = C(1,1,1,2)*2
 
  // mat[12] = C(2,2,0,0)
  // mat[13] = C(2,2,1,1)
  // mat[14] = C(2,2,2,2)
  // mat[15] = C(2,2,0,1)*2
  // mat[16] = C(2,2,0,2)*2
  // mat[17] = C(2,2,1,2)*2
 
  // mat[18] = C(0,1,0,0)
  // mat[19] = C(0,1,1,1)
  // mat[20] = C(0,1,2,2)
  // mat[21] = C(0,1,0,1)*2
  // mat[22] = C(0,1,0,2)*2
  // mat[23] = C(0,1,1,2)*2
 
  // mat[24] = C(0,2,0,0)
  // mat[25] = C(0,2,1,1)
  // mat[26] = C(0,2,2,2)
  // mat[27] = C(0,2,0,1)*2
  // mat[28] = C(0,2,0,2)*2
  // mat[29] = C(0,2,1,2)*2
 
  // mat[30] = C(1,2,0,0)
  // mat[31] = C(1,2,1,1)
  // mat[32] = C(1,2,2,2)
  // mat[33] = C(1,2,0,1)*2
  // mat[34] = C(1,2,0,2)*2
  // mat[35] = C(1,2,1,2)*2
 
  unsigned int index = 0;
  for (auto i : make_range(N))
    for (auto j : make_range(N))
      for (auto k : make_range(N))
        for (auto l : make_range(N))
        {
          if (i == j)
            mat[k == l ? i * ntens + k : i * ntens + k + nskip + l] += _vals[index];
          else
            // i!=j
            mat[k == l ? (nskip + i + j) * ntens + k : (nskip + i + j) * ntens + k + nskip + l] +=
                _vals[index]; // note the +=, which results in double-counting and is rectified
                              // below
          index++;
        }
 
  for (unsigned int i = 3; i < ntens; ++i)
    for (auto j : make_range(ntens))
      mat[i * ntens + j] /= 2.0; // because of double-counting above
 
  // use LAPACK to find the inverse
  MatrixTools::inverse(mat, ntens);
 
  // build the resulting rank-four tensor
  // using the inverse of the above algorithm
  index = 0;
  for (auto i : make_range(N))
    for (auto j : make_range(N))
      for (auto k : make_range(N))
        for (auto l : make_range(N))
        {
          if (i == j)
            result._vals[index] =
                k == l ? mat[i * ntens + k] : mat[i * ntens + k + nskip + l] / 2.0;
          else
            // i!=j
            result._vals[index] = k == l ? mat[(nskip + i + j) * ntens + k]
                                         : mat[(nskip + i + j) * ntens + k + nskip + l] / 2.0;
          index++;
        }
 
  return result;
}

◆ isIsotropic()

template<typename T >

bool RankFourTensorTempl< T >::isIsotropic ( ) const

checks if the tensor is isotropic

Definition at line 1130 of file RankFourTensorImplementation.h.

{
  // prerequisite is symmetry
  if (!isSymmetric())
    return false;
 
  // inspect shear components
  const T & mu = (*this)(0, 1, 0, 1);
  // ...diagonal
  if ((*this)(1, 2, 1, 2) != mu || (*this)(2, 0, 2, 0) != mu)
    return false;
  // ...off-diagonal
  if ((*this)(2, 0, 1, 2) != 0.0 || (*this)(0, 1, 1, 2) != 0.0 || (*this)(0, 1, 2, 0) != 0.0)
    return false;
 
  // off diagonal blocks in Voigt
  for (auto i : make_range(N))
    for (auto j : make_range(N))
      if (_vals[i * (N3 + N2) + ((j + 1) % N) * N + (j + 2) % N] != 0.0)
        return false;
 
  // top left block
  const T & K1 = (*this)(0, 0, 0, 0);
  const T & K2 = (*this)(0, 0, 1, 1);
  if (!MooseUtils::relativeFuzzyEqual(K1 - 4.0 * mu / 3.0, K2 + 2.0 * mu / 3.0))
    return false;
  if ((*this)(1, 1, 1, 1) != K1 || (*this)(2, 2, 2, 2) != K1)
    return false;
  for (auto i : make_range(1u, N))
    for (auto j : make_range(i))
      if ((*this)(i, i, j, j) != K2)
        return false;
 
  return true;
}

◆ isSymmetric()

template<typename T >

bool RankFourTensorTempl< T >::isSymmetric ( ) const

checks if the tensor is symmetric

Definition at line 1109 of file RankFourTensorImplementation.h.

{
  for (auto i : make_range(1u, N))
    for (auto j : make_range(i))
      for (auto k : make_range(1u, N))
        for (auto l : make_range(k))
        {
          // minor symmetries
          if ((*this)(i, j, k, l) != (*this)(j, i, k, l) ||
              (*this)(i, j, k, l) != (*this)(i, j, l, k))
            return false;
 
          // major symmetry
          if ((*this)(i, j, k, l) != (*this)(k, l, i, j))
            return false;
        }
  return true;
}

◆ L2norm()

template<typename T >

T RankFourTensorTempl< T >::L2norm ( ) const

sqrt(C_ijkl*C_ijkl)

Definition at line 268 of file RankFourTensorImplementation.h.

{
  T l2 = 0;
 
  for (auto i : make_range(N4))
    l2 += Utility::pow<2>(_vals[i]);
 
  using std::sqrt;
  return sqrt(l2);
}

◆ operator()() [1/2]

template<typename T >

T & RankFourTensorTempl< T >::operator()	(	unsigned int	i,
		unsigned int	j,
		unsigned int	k,
		unsigned int	l
	)

inline

Gets the value for the indices specified. Takes indices ranging from 0-2 for i, j, k, and l.

Definition at line 180 of file RankFourTensor.h.

  {
    return _vals[i * N3 + j * N2 + k * N + l];
  }

◆ operator()() [2/2]

template<typename T >

const T & RankFourTensorTempl< T >::operator()	(	unsigned int	i,
		unsigned int	j,
		unsigned int	k,
		unsigned int	l
	)		const

inline

Gets the value for the indices specified.

Takes indices ranging from 0-2 for i, j, k, and l. used for const

Definition at line 189 of file RankFourTensor.h.

  {
    return _vals[i * N3 + j * N2 + k * N + l];
  }

◆ operator*() [1/3]

template<typename T >

template<typename T2 >

auto RankFourTensorTempl< T >::operator* ( const RankFourTensorTempl< T2 > & a ) const -> RankFourTensorTempl<decltype(T() * T2())>

C_ijpq*a_pqkl.

Definition at line 234 of file RankFourTensorImplementation.h.

{
  typedef decltype(T() * T2()) ValueType;
  RankFourTensorTempl<ValueType> result;
 
  if constexpr (std::is_same_v<T, Real> && std::is_same_v<T2, Real>)
  {
    // C_ijkl = A_ijpq B_pqkl. With the row-major flat storage (linear index ij*N2 + kl, where
    // ij = i*N + j and kl = k*N + l) this is exactly an (N2 x N2) matrix product over the
    // flattened tensors. Dispatch to Eigen's vectorized kernel instead of the six-deep loop,
    // which recomputes the flat index via operator() ~3x per innermost iteration.
    const Eigen::Map<const Eigen::Matrix<Real, N2, N2, Eigen::RowMajor>> a(_vals);
    const Eigen::Map<const Eigen::Matrix<Real, N2, N2, Eigen::RowMajor>> bmat(b._vals);
    Eigen::Map<Eigen::Matrix<Real, N2, N2, Eigen::RowMajor>> r(result._vals);
    r.noalias() = a * bmat;
    return result;
  }
  else
  {
    for (auto i : make_range(N))
      for (auto j : make_range(N))
        for (auto k : make_range(N))
          for (auto l : make_range(N))
            for (auto p : make_range(N))
              for (auto q : make_range(N))
                result(i, j, k, l) += (*this)(i, j, p, q) * b(p, q, k, l);
 
    return result;
  }
}

◆ operator*() [2/3]

template<typename T >

template<typename T2 >

auto RankFourTensorTempl< T >::operator* ( const T2 & a ) const -> typename std::enable_if<libMesh::ScalarTraits<T2>::value, RankFourTensorTempl<decltype(T() * T2())>>::type

C_ijkl*a.

Definition at line 579 of file RankFourTensor.h.

{
  typedef decltype(T() * T2()) ValueType;
  RankFourTensorTempl<ValueType> result;
 
  for (auto i : libMesh::make_range(N4))
    result._vals[i] = _vals[i] * b;
 
  return result;
}

◆ operator*() [3/3]

template<typename T >

template<template< typename > class Tensor, typename T2 >

auto RankFourTensorTempl< T >::operator* ( const Tensor< T2 > & a ) const -> typename std::enable_if<TwoTensorMultTraits<Tensor, T2>::value, RankTwoTensorTempl<decltype(T() * T2())>>::type

C_ijkl*a_kl.

Definition at line 126 of file RankFourTensorImplementation.h.

{
  typedef decltype(T() * T2()) ValueType;
  RankTwoTensorTempl<ValueType> result;
 
  if constexpr (std::is_same_v<T, Real> && std::is_same_v<T2, Real>)
  {
    // result_ij = A_ijkl b_kl: an (N2 x N2) . (N2) matvec over the flat storage. Gather b into a
    // contiguous vector once (b's storage type is generic here), then use Eigen's kernel -- the
    // scalar loop below re-derives b's index with a div/mod on every inner iteration.
    Eigen::Matrix<Real, N2, 1> bvec;
    for (unsigned int kl = 0; kl < N2; ++kl)
      bvec(kl) = b(kl / N, kl % N);
    const Eigen::Map<const Eigen::Matrix<Real, N2, N2, Eigen::RowMajor>> a(_vals);
    Eigen::Map<Eigen::Matrix<Real, N2, 1>> r(result._coords);
    r.noalias() = a * bvec;
    return result;
  }
  else
  {
    unsigned int index = 0;
    for (unsigned int ij = 0; ij < N2; ++ij)
    {
      ValueType tmp = 0;
      for (unsigned int kl = 0; kl < N2; ++kl)
        tmp += _vals[index++] * b(kl / LIBMESH_DIM, kl % LIBMESH_DIM);
      result._coords[ij] = tmp;
    }
 
    return result;
  }
}

◆ operator*=()

template<typename T >

RankFourTensorTempl< T > & RankFourTensorTempl< T >::operator*= ( const T & a )

C_ijkl *= a.

Definition at line 163 of file RankFourTensorImplementation.h.

{
  for (auto i : make_range(N4))
    _vals[i] *= a;
  return *this;
}

◆ operator+()

template<typename T >

template<typename T2 >

auto RankFourTensorTempl< T >::operator+ ( const RankFourTensorTempl< T2 > & a ) const -> RankFourTensorTempl<decltype(T() + T2())>

C_ijkl + a_ijkl.

Definition at line 191 of file RankFourTensorImplementation.h.

{
  RankFourTensorTempl<decltype(T() + T2())> result;
  for (auto i : make_range(N4))
    result._vals[i] = _vals[i] + b._vals[i];
  return result;
}

◆ operator+=()

template<typename T >

RankFourTensorTempl< T > & RankFourTensorTempl< T >::operator+= ( const RankFourTensorTempl< T > & a )

C_ijkl += a_ijkl for all i, j, k, l.

Definition at line 181 of file RankFourTensorImplementation.h.

{
  for (auto i : make_range(N4))
    _vals[i] += a._vals[i];
  return *this;
}

◆ operator-() [1/2]

template<typename T >

RankFourTensorTempl< T > RankFourTensorTempl< T >::operator- ( ) const

-C_ijkl

Definition at line 223 of file RankFourTensorImplementation.h.

{
  RankFourTensorTempl<T> result;
  for (auto i : make_range(N4))
    result._vals[i] = -_vals[i];
  return result;
}

◆ operator-() [2/2]

template<typename T >

template<typename T2 >

auto RankFourTensorTempl< T >::operator- ( const RankFourTensorTempl< T2 > & a ) const -> RankFourTensorTempl<decltype(T() - T2())>

C_ijkl - a_ijkl.

Definition at line 212 of file RankFourTensorImplementation.h.

{
  RankFourTensorTempl<decltype(T() - T2())> result;
  for (auto i : make_range(N4))
    result._vals[i] = _vals[i] - b._vals[i];
  return result;
}

◆ operator-=()

template<typename T >

RankFourTensorTempl< T > & RankFourTensorTempl< T >::operator-= ( const RankFourTensorTempl< T > & a )

C_ijkl -= a_ijkl.

Definition at line 202 of file RankFourTensorImplementation.h.

{
  for (auto i : make_range(N4))
    _vals[i] -= a._vals[i];
  return *this;
}

◆ operator/()

template<typename T >

template<typename T2 >

auto RankFourTensorTempl< T >::operator/ ( const T2 & a ) const -> typename std::enable_if<libMesh::ScalarTraits<T2>::value, RankFourTensorTempl<decltype(T() / T2())>>::type

C_ijkl/a.

Definition at line 595 of file RankFourTensor.h.

{
  RankFourTensorTempl<decltype(T() / T2())> result;
  for (auto i : libMesh::make_range(N4))
    result._vals[i] = _vals[i] / b;
  return result;
}

◆ operator/=()

template<typename T >

RankFourTensorTempl< T > & RankFourTensorTempl< T >::operator/= ( const T & a )

C_ijkl /= a for all i, j, k, l.

Definition at line 172 of file RankFourTensorImplementation.h.

{
  for (auto i : make_range(N4))
    _vals[i] /= a;
  return *this;
}

◆ operator=() [1/2]

template<typename T >

RankFourTensorTempl< T > & RankFourTensorTempl< T >::operator= ( const RankFourTensorTempl< T > & a )

default

copies values from a into this tensor

◆ operator=() [2/2]

template<typename T >

template<typename Scalar >

std::enable_if< libMesh::ScalarTraits< Scalar >::value, RankFourTensorTempl & >::type RankFourTensorTempl< T >::operator= ( const Scalar & libmesh_dbg_varp )

inline

Assignment-from-scalar operator.

Used only to zero out the tensor.

Returns: A reference to *this.

Definition at line 219 of file RankFourTensor.h.

  {
    libmesh_assert_equal_to(p, Scalar(0));
    this->zero();
    return *this;
  }

◆ print()

template<typename T >

void RankFourTensorTempl< T >::print ( std::ostream & stm = Moose::out ) const

Print the rank four tensor.

Definition at line 444 of file RankFourTensorImplementation.h.

{
  for (auto i : make_range(N))
    for (auto j : make_range(N))
    {
      stm << "i = " << i << " j = " << j << '\n';
      for (auto k : make_range(N))
      {
        for (auto l : make_range(N))
          stm << std::setw(15) << (*this)(i, j, k, l) << " ";
 
        stm << '\n';
      }
    }
 
  stm << std::flush;
}

◆ printReal()

template<typename T >

void RankFourTensorTempl< T >::printReal ( std::ostream & stm = Moose::out ) const

Print the values of the rank four tensor.

Definition at line 464 of file RankFourTensorImplementation.h.

{
  for (unsigned int i = 0; i < N; ++i)
    for (unsigned int j = 0; j < N; ++j)
    {
      stm << "i = " << i << " j = " << j << '\n';
      for (unsigned int k = 0; k < N; ++k)
      {
        for (unsigned int l = 0; l < N; ++l)
          stm << std::setw(15) << MetaPhysicL::raw_value((*this)(i, j, k, l)) << " ";
 
        stm << '\n';
      }
    }
 
  stm << std::flush;
}

◆ rotate()

template<typename T >

void RankFourTensorTempl< T >::rotate ( const TypeTensor< T > & R )

Rotate the tensor using C_ijkl = R_im R_jn R_ko R_lp C_mnop.

Definition at line 412 of file RankFourTensorImplementation.h.

{
  RankFourTensorTempl<T> old = *this;
 
  unsigned int index = 0;
  for (auto i : make_range(N))
    for (auto j : make_range(N))
      for (auto k : make_range(N))
        for (auto l : make_range(N))
        {
          unsigned int index2 = 0;
          T sum = 0.0;
          for (auto m : make_range(N))
          {
            const T & a = R(i, m);
            for (auto n : make_range(N))
            {
              const T & ab = a * R(j, n);
              for (auto o : make_range(N))
              {
                const T & abc = ab * R(k, o);
                for (auto p : make_range(N))
                  sum += abc * R(l, p) * old._vals[index2++];
              }
            }
          }
          _vals[index++] = sum;
        }
}

◆ singleProductI()

template<typename T >

RankFourTensorTempl< T > RankFourTensorTempl< T >::singleProductI ( const RankTwoTensorTempl< T > & A ) const

Calculates C_mjkl A_im.

Definition at line 1045 of file RankFourTensorImplementation.h.

{
  RankFourTensorTempl<T> R;
 
  for (unsigned int i = 0; i < N; i++)
    for (unsigned int j = 0; j < N; j++)
      for (unsigned int k = 0; k < N; k++)
        for (unsigned int l = 0; l < N; l++)
          for (unsigned int m = 0; m < N; m++)
            R(i, j, k, l) += (*this)(m, j, k, l) * A(i, m);
 
  return R;
}

◆ singleProductJ()

template<typename T >

RankFourTensorTempl< T > RankFourTensorTempl< T >::singleProductJ ( const RankTwoTensorTempl< T > & A ) const

Calculates C_imkl A_jm.

Definition at line 1061 of file RankFourTensorImplementation.h.

{
  RankFourTensorTempl<T> R;
 
  for (unsigned int i = 0; i < N; i++)
    for (unsigned int j = 0; j < N; j++)
      for (unsigned int k = 0; k < N; k++)
        for (unsigned int l = 0; l < N; l++)
          for (unsigned int m = 0; m < N; m++)
            R(i, j, k, l) += (*this)(i, m, k, l) * A(j, m);
 
  return R;
}

◆ singleProductK()

template<typename T >

RankFourTensorTempl< T > RankFourTensorTempl< T >::singleProductK ( const RankTwoTensorTempl< T > & A ) const

Calculates C_ijml A_km.

Definition at line 1077 of file RankFourTensorImplementation.h.

{
  RankFourTensorTempl<T> R;
 
  for (unsigned int i = 0; i < N; i++)
    for (unsigned int j = 0; j < N; j++)
      for (unsigned int k = 0; k < N; k++)
        for (unsigned int l = 0; l < N; l++)
          for (unsigned int m = 0; m < N; m++)
            R(i, j, k, l) += (*this)(i, j, m, l) * A(k, m);
 
  return R;
}

◆ singleProductL()

template<typename T >

RankFourTensorTempl< T > RankFourTensorTempl< T >::singleProductL ( const RankTwoTensorTempl< T > & A ) const

Calculates C_ijkm A_lm.

Definition at line 1093 of file RankFourTensorImplementation.h.

{
  RankFourTensorTempl<T> R;
 
  for (unsigned int i = 0; i < N; i++)
    for (unsigned int j = 0; j < N; j++)
      for (unsigned int k = 0; k < N; k++)
        for (unsigned int l = 0; l < N; l++)
          for (unsigned int m = 0; m < N; m++)
            R(i, j, k, l) += (*this)(i, j, k, m) * A(l, m);
 
  return R;
}

◆ sum3x1()

template<typename T >

libMesh::VectorValue< T > RankFourTensorTempl< T >::sum3x1 ( ) const

Calculates the vector a[i] = sum over j Ciijj for i and j varying from 0 to 2.

Definition at line 957 of file RankFourTensorImplementation.h.

{
  // used for volumetric locking correction
  libMesh::VectorValue<T> a(3);
  a(0) = (*this)(0, 0, 0, 0) + (*this)(0, 0, 1, 1) + (*this)(0, 0, 2, 2); // C0000 + C0011 + C0022
  a(1) = (*this)(1, 1, 0, 0) + (*this)(1, 1, 1, 1) + (*this)(1, 1, 2, 2); // C1100 + C1111 + C1122
  a(2) = (*this)(2, 2, 0, 0) + (*this)(2, 2, 1, 1) + (*this)(2, 2, 2, 2); // C2200 + C2211 + C2222
  return a;
}

◆ sum3x3()

template<typename T >

T RankFourTensorTempl< T >::sum3x3 ( ) const

Calculates the sum of Ciijj for i and j varying from 0 to 2.

Definition at line 945 of file RankFourTensorImplementation.h.

{
  // used in the volumetric locking correction
  T sum = 0;
  for (auto i : make_range(N))
    for (auto j : make_range(N))
      sum += (*this)(i, i, j, j);
  return sum;
}

◆ surfaceFillFromInputVector()

template<typename T >

void RankFourTensorTempl< T >::surfaceFillFromInputVector ( const std::vector< T > & input )

Fills the tensor entries ignoring the last dimension (ie, C_ijkl=0 if any of i, j, k, or l = 3).

Fill method depends on size of input Input size = 2. Then C_1111 = C_2222 = input[0], and C_1122 = input[1], and C_1212 = (input[0]

input[1])/2, and C_ijkl = C_jikl = C_ijlk = C_klij, and C_1211 = C_1222 = 0. Input size = 9. Then C_1111 = input[0], C_1112 = input[1], C_1122 = input[3], C_1212 = input[4], C_1222 = input[5], C_1211 = input[6] C_2211 = input[7], C_2212 = input[8], C_2222 = input[9] and C_ijkl = C_jikl = C_ijlk

Definition at line 530 of file RankFourTensorImplementation.h.

{
  zero();
 
  if (input.size() == 9)
  {
    // then fill from vector C_1111, C_1112, C_1122, C_1212, C_1222, C_1211, C_2211, C_2212, C_2222
    (*this)(0, 0, 0, 0) = input[0];
    (*this)(0, 0, 0, 1) = input[1];
    (*this)(0, 0, 1, 1) = input[2];
    (*this)(0, 1, 0, 1) = input[3];
    (*this)(0, 1, 1, 1) = input[4];
    (*this)(0, 1, 0, 0) = input[5];
    (*this)(1, 1, 0, 0) = input[6];
    (*this)(1, 1, 0, 1) = input[7];
    (*this)(1, 1, 1, 1) = input[8];
 
    // fill in remainders from C_ijkl = C_ijlk = C_jikl
    (*this)(0, 0, 1, 0) = (*this)(0, 0, 0, 1);
    (*this)(0, 1, 1, 0) = (*this)(0, 1, 0, 1);
    (*this)(1, 0, 0, 0) = (*this)(0, 1, 0, 0);
    (*this)(1, 0, 0, 1) = (*this)(0, 1, 0, 1);
    (*this)(1, 0, 1, 1) = (*this)(0, 1, 1, 1);
    (*this)(1, 0, 0, 0) = (*this)(0, 1, 0, 0);
    (*this)(1, 1, 1, 0) = (*this)(1, 1, 0, 1);
  }
  else if (input.size() == 2)
  {
    // only two independent constants, C_1111 and C_1122
    (*this)(0, 0, 0, 0) = input[0];
    (*this)(0, 0, 1, 1) = input[1];
    // use symmetries
    (*this)(1, 1, 1, 1) = (*this)(0, 0, 0, 0);
    (*this)(1, 1, 0, 0) = (*this)(0, 0, 1, 1);
    (*this)(0, 1, 0, 1) = 0.5 * ((*this)(0, 0, 0, 0) - (*this)(0, 0, 1, 1));
    (*this)(1, 0, 0, 1) = (*this)(0, 1, 0, 1);
    (*this)(0, 1, 1, 0) = (*this)(0, 1, 0, 1);
    (*this)(1, 0, 1, 0) = (*this)(0, 1, 0, 1);
  }
  else
    mooseError("Please provide correct number of inputs for surface RankFourTensorTempl<T> "
               "initialization.");
}

◆ transposeIj()

template<typename T >

RankFourTensorTempl< T > RankFourTensorTempl< T >::transposeIj ( ) const

Transpose the tensor by swapping the first two indeces.

Returns: C_jikl

Definition at line 500 of file RankFourTensorImplementation.h.

{
  RankFourTensorTempl<T> result;
 
  for (auto i : make_range(N))
    for (auto j : make_range(N))
      for (auto k : make_range(N))
        for (auto l : make_range(N))
          result(i, j, k, l) = (*this)(j, i, k, l);
 
  return result;
}

◆ transposeKl()

template<typename T >

RankFourTensorTempl< T > RankFourTensorTempl< T >::transposeKl ( ) const

Transpose the tensor by swapping the last two indeces.

Returns: C_ijlk

Definition at line 515 of file RankFourTensorImplementation.h.

{
  RankFourTensorTempl<T> result;
 
  for (auto i : make_range(N))
    for (auto j : make_range(N))
      for (auto k : make_range(N))
        for (auto l : make_range(N))
          result(i, j, k, l) = (*this)(i, j, l, k);
 
  return result;
}

◆ transposeMajor()

template<typename T >

RankFourTensorTempl< T > RankFourTensorTempl< T >::transposeMajor ( ) const

Transpose the tensor by swapping the first pair with the second pair of indices.

Returns: C_klji

Definition at line 484 of file RankFourTensorImplementation.h.

{
  RankFourTensorTempl<T> result;
 
  unsigned int index = 0;
  for (auto i : make_range(N))
    for (auto j : make_range(N))
      for (auto k : make_range(N))
        for (auto l : make_range(N))
          result._vals[index++] = _vals[k * N3 + i * N + j + l * N2];
 
  return result;
}

◆ tripleProductIjk()

template<typename T >

RankFourTensorTempl< T > RankFourTensorTempl< T >::tripleProductIjk	(	const RankTwoTensorTempl< T > &	A,
		const RankTwoTensorTempl< T > &	B,
		const RankTwoTensorTempl< T > &	C
	)		const

Calculates C_mntl A_im B_jn C_kt.

Definition at line 1026 of file RankFourTensorImplementation.h.

{
  RankFourTensorTempl<T> R;
  for (unsigned int i = 0; i < N; i++)
    for (unsigned int j = 0; j < N; j++)
      for (unsigned int k = 0; k < N; k++)
        for (unsigned int l = 0; l < N; l++)
          for (unsigned int m = 0; m < N; m++)
            for (unsigned int n = 0; n < N; n++)
              for (unsigned int t = 0; t < N; t++)
                R(i, j, k, l) += (*this)(m, n, t, l) * A(i, m) * B(j, n) * C(k, t);
 
  return R;
}

◆ tripleProductIjl()

template<typename T >

RankFourTensorTempl< T > RankFourTensorTempl< T >::tripleProductIjl	(	const RankTwoTensorTempl< T > &	A,
		const RankTwoTensorTempl< T > &	B,
		const RankTwoTensorTempl< T > &	C
	)		const

Calculates C_mnkt A_im B_jn C_lt.

Definition at line 1007 of file RankFourTensorImplementation.h.

{
  RankFourTensorTempl<T> R;
  for (unsigned int i = 0; i < N; i++)
    for (unsigned int j = 0; j < N; j++)
      for (unsigned int k = 0; k < N; k++)
        for (unsigned int l = 0; l < N; l++)
          for (unsigned int m = 0; m < N; m++)
            for (unsigned int n = 0; n < N; n++)
              for (unsigned int t = 0; t < N; t++)
                R(i, j, k, l) += (*this)(m, n, k, t) * A(i, m) * B(j, n) * C(l, t);
 
  return R;
}

◆ tripleProductIkl()

template<typename T >

RankFourTensorTempl< T > RankFourTensorTempl< T >::tripleProductIkl	(	const RankTwoTensorTempl< T > &	A,
		const RankTwoTensorTempl< T > &	B,
		const RankTwoTensorTempl< T > &	C
	)		const

Calculates C_mjnt A_im B_kn C_lt.

Definition at line 988 of file RankFourTensorImplementation.h.

{
  RankFourTensorTempl<T> R;
  for (unsigned int i = 0; i < N; i++)
    for (unsigned int j = 0; j < N; j++)
      for (unsigned int k = 0; k < N; k++)
        for (unsigned int l = 0; l < N; l++)
          for (unsigned int m = 0; m < N; m++)
            for (unsigned int n = 0; n < N; n++)
              for (unsigned int t = 0; t < N; t++)
                R(i, j, k, l) += (*this)(m, j, n, t) * A(i, m) * B(k, n) * C(l, t);
 
  return R;
}

◆ tripleProductJkl()

template<typename T >

RankFourTensorTempl< T > RankFourTensorTempl< T >::tripleProductJkl	(	const RankTwoTensorTempl< T > &	A,
		const RankTwoTensorTempl< T > &	B,
		const RankTwoTensorTempl< T > &	C
	)		const

Calculates C_imnt A_jm B_kn C_lt.

Definition at line 969 of file RankFourTensorImplementation.h.

{
  RankFourTensorTempl<T> R;
  for (unsigned int i = 0; i < N; i++)
    for (unsigned int j = 0; j < N; j++)
      for (unsigned int k = 0; k < N; k++)
        for (unsigned int l = 0; l < N; l++)
          for (unsigned int m = 0; m < N; m++)
            for (unsigned int n = 0; n < N; n++)
              for (unsigned int t = 0; t < N; t++)
                R(i, j, k, l) += (*this)(i, m, n, t) * A(j, m) * B(k, n) * C(l, t);
 
  return R;
}

◆ zero()

template<typename T >

void RankFourTensorTempl< T >::zero ( )

Zeros out the tensor.

Definition at line 117 of file RankFourTensorImplementation.h.

{
  for (auto i : make_range(N4))
    _vals[i] = 0.0;
}

Referenced by MathUtils::mooseSetToZero< ADRankFourTensor >(), MathUtils::mooseSetToZero< RankFourTensor >(), and RankFourTensorTempl< T >::operator=().

Friends And Related Symbol Documentation

◆ dataLoad

template<typename T >

template<class T2 >

void dataLoad	(	std::istream &	,
		RankFourTensorTempl< T2 > &	,
		void *
	)

friend

◆ dataStore

template<typename T >

template<class T2 >

void dataStore	(	std::ostream &	,
		RankFourTensorTempl< T2 > &	,
		void *
	)

friend

◆ operator<<

template<typename T >

std::ostream & operator<<	(	std::ostream &	os,
		const RankFourTensorTempl< T > &	t
	)

friend

Definition at line 200 of file RankFourTensor.h.

  {
    t.print(os);
    return os;
  }

◆ RankFourTensorTempl

template<typename T >

template<typename T2 >

friend class RankFourTensorTempl

friend

Definition at line 532 of file RankFourTensor.h.

◆ RankThreeTensorTempl

template<typename T >

template<typename T2 >

friend class RankThreeTensorTempl

friend

Definition at line 534 of file RankFourTensor.h.

◆ RankTwoTensorTempl

template<typename T >

template<typename T2 >

friend class RankTwoTensorTempl

friend

Definition at line 530 of file RankFourTensor.h.

Member Data Documentation

◆ _vals

template<typename T >

T RankFourTensorTempl< T >::_vals[N4]

protected

The values of the rank-four tensor stored by index=(((i * LIBMESH_DIM + j) * LIBMESH_DIM + k) * LIBMESH_DIM + l)

Definition at line 435 of file RankFourTensor.h.

◆ N

template<typename T >

constexpr unsigned int RankFourTensorTempl< T >::N = Moose::dim

staticconstexpr

tensor dimension and powers of the dimension

Definition at line 69 of file RankFourTensor.h.

Referenced by RankFourTensorTempl< T >::operator()(), and RankFourTensorTempl< T >::operator()().

◆ N2

template<typename T >

constexpr unsigned int RankFourTensorTempl< T >::N2 = N * N

staticconstexpr

Definition at line 70 of file RankFourTensor.h.

Referenced by RankFourTensorTempl< T >::operator()(), and RankFourTensorTempl< T >::operator()().

◆ N3

template<typename T >

constexpr unsigned int RankFourTensorTempl< T >::N3 = N * N * N

staticconstexpr

Definition at line 71 of file RankFourTensor.h.

Referenced by RankFourTensorTempl< T >::operator()(), and RankFourTensorTempl< T >::operator()().

◆ N4

template<typename T >

constexpr unsigned int RankFourTensorTempl< T >::N4 = N * N * N * N

staticconstexpr

Definition at line 72 of file RankFourTensor.h.

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

include/utils/RankFourTensor.h
include/utils/RankFourTensorImplementation.h

Classes

Public Types

Public Member Functions

Static Public Member Functions

Static Public Attributes

Protected Member Functions

Protected Attributes

Friends

Detailed Description

Member Typedef Documentation

◆ index_type

◆ value_type

Member Enumeration Documentation

◆ FillMethod

◆ InitMethod

Constructor & Destructor Documentation

◆ RankFourTensorTempl() [1/6]

◆ RankFourTensorTempl() [2/6]

◆ RankFourTensorTempl() [3/6]

◆ RankFourTensorTempl() [4/6]

◆ RankFourTensorTempl() [5/6]

◆ RankFourTensorTempl() [6/6]

Member Function Documentation

◆ contraction()

◆ contractionIj()

◆ contractionKl()

◆ fillAntisymmetricFromInputVector()

◆ fillAntisymmetricIsotropic()

◆ fillAntisymmetricIsotropicFromInputVector()

◆ fillAxisymmetricRZFromInputVector()

◆ fillFromInputVector()

◆ fillGeneralFromInputVector()

◆ fillGeneralIsotropic()

◆ fillGeneralIsotropicFromInputVector()

◆ fillGeneralOrthotropicFromInputVector()

◆ fillMethodEnum()

◆ fillPrincipalFromInputVector()

◆ fillSymmetric21FromInputVector()

◆ fillSymmetric9FromInputVector()

◆ fillSymmetricIsotropic()

◆ fillSymmetricIsotropicEandNu()

◆ fillSymmetricIsotropicEandNuFromInputVector()

◆ fillSymmetricIsotropicFromInputVector()

◆ Identity()

◆ IdentityDeviatoric()

◆ IdentityFour()

◆ innerProductTranspose()

◆ inverse()

◆ invSymm() [1/2]

◆ invSymm() [2/2]

◆ isIsotropic()

◆ isSymmetric()

◆ L2norm()

◆ operator()() [1/2]

◆ operator()() [2/2]

◆ operator*() [1/3]

◆ operator*() [2/3]

◆ operator*() [3/3]

◆ operator*=()

◆ operator+()

◆ operator+=()

◆ operator-() [1/2]

◆ operator-() [2/2]

◆ operator-=()

◆ operator/()

◆ operator/=()

◆ operator=() [1/2]

◆ operator=() [2/2]

◆ print()

◆ printReal()

◆ rotate()

◆ singleProductI()

◆ singleProductJ()

◆ singleProductK()

◆ singleProductL()

◆ sum3x1()

◆ sum3x3()

◆ surfaceFillFromInputVector()

◆ transposeIj()

◆ transposeKl()