v0.8.13
Classes | Typedefs | Enumerations | Functions | Variables
EshelbianPlasticity Namespace Reference

Classes

struct  BcDisp
 
struct  BcLoadScale
 
struct  DataAtIntegrationPts
 
struct  EshelbianCore
 
struct  FaceRule
 
struct  FeTractionFreeBc
 
struct  HMHStVenantKirchhoff
 
struct  OpAssembleBasic
 
struct  OpAssembleFace
 
struct  OpAssembleVolume
 
struct  OpCalculateRotationAndSpatialGradient
 
struct  OpDispBc
 
struct  OpHMHHStVenantKirchhoff
 
struct  OpJacobian
 
struct  OpPostProcDataStructure
 
struct  OpSpatialConsistency
 
struct  OpSpatialConsistency_dP_domega
 
struct  OpSpatialEquilibrium
 
struct  OpSpatialEquilibrium_dw_dP
 
struct  OpSpatialPhysical
 
struct  OpSpatialPhysical_du_domega
 
struct  OpSpatialPhysical_du_dP
 
struct  OpSpatialPhysical_du_du
 
struct  OpSpatialRotation
 
struct  OpSpatialRotation_domega_domega
 
struct  OpSpatialRotation_domega_dP
 
struct  PhysicalEquations
 
struct  VolRule
 Set integration rule on element. More...
 

Typedefs

typedef boost::shared_ptr< MatrixDouble > MatrixPtr
 
typedef boost::shared_ptr< VectorDouble > VectorPtr
 
using EntData = DataForcesAndSourcesCore::EntData
 
using UserDataOperator = ForcesAndSourcesCore::UserDataOperator
 
using VolUserDataOperator = VolumeElementForcesAndSourcesCore::UserDataOperator
 
using FaceUserDataOperator = FaceElementForcesAndSourcesCore::UserDataOperator
 
typedef std::vector< BcDispBcDispVec
 
typedef std::vector< Range > TractionFreeBc
 

Enumerations

enum  EshelbianInterfaces { CORE_ESHELBIAN_INTERFACE = 1 << 0 }
 

Functions

FTensor::Tensor3< FTensor::PackPtr< double *, 1 >, 3, 3, 3 > getFTensor3FromMat (MatrixDouble &m)
 
template<typename T >
FTensor::Tensor2< double, 3, 3 > getRotationFormVector (FTensor::Tensor1< T, 3 > &t_omega)
 
template<typename T >
FTensor::Tensor3< double, 3, 3, 3 > getDiffRotationFormVector (FTensor::Tensor1< T, 3 > &t_omega)
 

Variables

static const MOFEMuuid IDD_MOFEMEshelbianCrackInterface
 

Typedef Documentation

◆ BcDispVec

typedef std::vector<BcDisp> EshelbianPlasticity::BcDispVec

Definition at line 334 of file EshelbianPlasticity.hpp.

◆ EntData

using EshelbianPlasticity::EntData = typedef DataForcesAndSourcesCore::EntData

◆ FaceUserDataOperator

using EshelbianPlasticity::FaceUserDataOperator = typedef FaceElementForcesAndSourcesCore::UserDataOperator

Definition at line 105 of file EshelbianPlasticity.hpp.

◆ MatrixPtr

typedef boost::shared_ptr<MatrixDouble> EshelbianPlasticity::MatrixPtr

Definition at line 16 of file EshelbianPlasticity.hpp.

◆ TractionFreeBc

typedef std::vector<Range> EshelbianPlasticity::TractionFreeBc

Definition at line 335 of file EshelbianPlasticity.hpp.

◆ UserDataOperator

using EshelbianPlasticity::UserDataOperator = typedef ForcesAndSourcesCore::UserDataOperator

◆ VectorPtr

typedef boost::shared_ptr<VectorDouble> EshelbianPlasticity::VectorPtr

Definition at line 17 of file EshelbianPlasticity.hpp.

◆ VolUserDataOperator

using EshelbianPlasticity::VolUserDataOperator = typedef VolumeElementForcesAndSourcesCore::UserDataOperator

Definition at line 104 of file EshelbianPlasticity.hpp.

Enumeration Type Documentation

◆ EshelbianInterfaces

Enumerator
CORE_ESHELBIAN_INTERFACE 

Definition at line 11 of file EshelbianPlasticity.hpp.

Function Documentation

◆ getDiffRotationFormVector()

template<typename T >
FTensor::Tensor3<double, 3, 3, 3> EshelbianPlasticity::getDiffRotationFormVector ( FTensor::Tensor1< T, 3 > &  t_omega)

Definition at line 59 of file EshelbianPlasticity.hpp.

59  {
60 
65 
66  const double angle = sqrt(t_omega(i) * t_omega(i));
67 
68  if (fabs(angle) == 0) {
70  auto t_R = getRotationFormVector(t_omega);
71  t_diff_R(i, j, k) = levi_civita(i, l, k) * t_R(l, j);
72  return t_diff_R;
73  }
74 
77  t_Omega(i, j) = levi_civita(i, j, k) * t_omega(k);
78 
79  const double angle2 = angle * angle;
80  const double ss = sin(angle);
81  const double a = ss / angle;
82  const double ss_2 = sin(angle / 2.);
83  const double b = 2 * ss_2 * ss_2 / angle2;
84  t_diff_R(i, j, k) = 0;
85  t_diff_R(i, j, k) = a * levi_civita(i, j, k);
86  t_diff_R(i, j, k) += b * (t_Omega(i, l) * levi_civita(l, j, k) +
87  levi_civita(i, l, k) * t_Omega(l, j));
88 
89  const double cc = cos(angle);
90  const double cc_2 = cos(angle / 2.);
91  const double diff_a = (angle * cc - ss) / angle2;
92 
93  const double diff_b =
94  (2. * angle * ss_2 * cc_2 - 4. * ss_2 * ss_2) / (angle2 * angle);
95  t_diff_R(i, j, k) += diff_a * t_Omega(i, j) * (t_omega(k) / angle);
96  t_diff_R(i, j, k) +=
97  diff_b * t_Omega(i, l) * t_Omega(l, j) * (t_omega(k) / angle);
98 
99  return t_diff_R;
100 }
FTensor::Tensor2< double, 3, 3 > getRotationFormVector(FTensor::Tensor1< T, 3 > &t_omega)
constexpr std::enable_if<(Dim0<=2 &&Dim1<=2), Tensor2_Expr< Levi_Civita< T >, T, Dim0, Dim1, i, j > >::type levi_civita(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
levi_civita functions to make for easy adhoc use

◆ getFTensor3FromMat()

FTensor::Tensor3<FTensor::PackPtr<double *, 1>, 3, 3, 3> EshelbianPlasticity::getFTensor3FromMat ( MatrixDouble &  m)

Definition at line 20 of file EshelbianPlasticity.hpp.

20  {
22  &m(0, 0), &m(1, 0), &m(2, 0), &m(3, 0), &m(4, 0), &m(5, 0), &m(6, 0),
23  &m(7, 0), &m(8, 0), &m(9, 0), &m(10, 0), &m(11, 0), &m(12, 0), &m(13, 0),
24  &m(14, 0), &m(15, 0), &m(16, 0), &m(17, 0), &m(18, 0), &m(19, 0),
25  &m(20, 0), &m(21, 0), &m(22, 0), &m(23, 0), &m(24, 0), &m(25, 0),
26  &m(26, 0));
27 }

◆ getRotationFormVector()

template<typename T >
FTensor::Tensor2<double, 3, 3> EshelbianPlasticity::getRotationFormVector ( FTensor::Tensor1< T, 3 > &  t_omega)

Definition at line 31 of file EshelbianPlasticity.hpp.

31  {
36  t_R(i, j) = 0;
37  t_R(0, 0) = 1;
38  t_R(1, 1) = 1;
39  t_R(2, 2) = 1;
40 
41  const double angle = sqrt(t_omega(i) * t_omega(i));
42  if (fabs(angle) == 0) {
43  return t_R;
44  }
45 
47  t_Omega(i, j) = levi_civita(i, j, k) * t_omega(k);
48  const double a = sin(angle) / angle;
49  const double ss_2 = sin(angle / 2.);
50  const double b = 2 * ss_2 * ss_2 / (angle * angle);
51  t_R(i, j) += a * t_Omega(i, j);
52  t_R(i, j) += b * t_Omega(i, k) * t_Omega(k, j);
53 
54  return t_R;
55 }
constexpr std::enable_if<(Dim0<=2 &&Dim1<=2), Tensor2_Expr< Levi_Civita< T >, T, Dim0, Dim1, i, j > >::type levi_civita(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
levi_civita functions to make for easy adhoc use

Variable Documentation

◆ IDD_MOFEMEshelbianCrackInterface

const MOFEMuuid EshelbianPlasticity::IDD_MOFEMEshelbianCrackInterface
static
Initial value:

Definition at line 13 of file EshelbianPlasticity.hpp.