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

Classes

struct  EpElementBase
 
struct  EpElement
 
struct  EpElement< BasicMethod >
 
struct  EpElement< FEMethod >
 
struct  EpFEMethod
 
struct  DataAtIntegrationPts
 
struct  PhysicalEquations
 
struct  BcDisp
 
struct  BcRot
 
struct  TractionBc
 
struct  OpJacobian
 
struct  OpAssembleBasic
 
struct  OpAssembleVolume
 
struct  OpAssembleFace
 
struct  OpCalculateRotationAndSpatialGradient
 
struct  OpSpatialEquilibrium
 
struct  OpSpatialRotation
 
struct  OpSpatialPhysical
 
struct  OpSpatialConsistencyP
 
struct  OpSpatialConsistencyBubble
 
struct  OpSpatialConsistencyDivTerm
 
struct  OpDispBc
 
struct  OpDispBc_dx
 
struct  OpRotationBc
 
struct  OpRotationBc_dx
 
struct  OpTractionBc
 
struct  FeTractionBc
 
struct  EpElement< FeTractionBc >
 
struct  OpSpatialEquilibrium_dw_dP
 
struct  OpSpatialEquilibrium_dw_dw
 
struct  OpSpatialPhysical_du_du
 
struct  OpSpatialPhysical_du_dP
 
struct  OpSpatialPhysical_du_dBubble
 
struct  OpSpatialPhysical_du_domega
 
struct  OpSpatialPhysical_du_dx
 
struct  OpSpatialRotation_domega_dP
 
struct  OpSpatialRotation_domega_dBubble
 
struct  OpSpatialRotation_domega_domega
 
struct  OpSpatialConsistency_dP_domega
 
struct  OpSpatialConsistency_dBubble_domega
 
struct  OpPostProcDataStructure
 
struct  OpSpatialPrj
 
struct  OpSpatialPrj_dx_dx
 
struct  OpSpatialPrj_dx_dw
 
struct  OpSpatialSchurBegin
 
struct  OpSpatialPreconditionMass
 
struct  OpSpatialSchurEnd
 
struct  OpCalculateStrainEnergy
 
struct  OpL2Transform
 
struct  EshelbianCore
 
struct  OpHMHH
 
struct  HMHStVenantKirchhoff
 
struct  HMHPMooneyRivlinWriggersEq63
 
struct  TensorTypeExtractor
 
struct  TensorTypeExtractor< FTensor::PackPtr< T, I > >
 
struct  VolRule
 Set integration rule on element. More...
 
struct  FaceRule
 
struct  CGGUserPolynomialBase
 

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< BcRotBcRotVec
 
typedef std::vector< Range > TractionFreeBc
 
typedef std::vector< TractionBcTractionBcVec
 

Enumerations

enum  EshelbianInterfaces { CORE_ESHELBIAN_INTERFACE = 1 << 0 }
 

Functions

MoFEMErrorCode CGG_BubbleBase_MBTET (const int p, const double *N, const double *diffN, FTensor::Tensor2< FTensor::PackPtr< double *, 9 >, 3, 3 > &phi, const int gdim)
 Calculate CGGT tonsorial bubble base. More...
 
FTensor::Tensor3< FTensor::PackPtr< double *, 1 >, 3, 3, 3 > getFTensor3FromMat (MatrixDouble &m)
 
double f (const double v)
 
double d_f (const double v)
 
double dd_f (const double v)
 
template<typename T >
static FTensor::Tensor2< typename TensorTypeExtractor< T >::Type, 3, 3 > get_rotation_form_vector (FTensor::Tensor1< T, 3 > &t_omega)
 
template<typename T >
static FTensor::Tensor3< typename TensorTypeExtractor< T >::Type, 3, 3, 3 > get_diff_rotation_form_vector (FTensor::Tensor1< T, 3 > &t_omega)
 
template<typename T1 , typename T2 , typename T3 >
MoFEMErrorCode get_eigen_val_and_proj_lapack (T1 &X, T2 &eig, T3 &eig_vec)
 
bool is_eq (const double &a, const double &b)
 
template<typename T >
int get_uniq_nb (T &ec)
 
template<typename T1 , typename T2 >
void sort_eigen_vals (T1 &eig, T2 &eigen_vec)
 

Variables

static const MOFEMuuid IDD_MOFEMEshelbianCrackInterface
 

Typedef Documentation

◆ BcDispVec

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

Definition at line 554 of file EshelbianPlasticity.hpp.

◆ BcRotVec

typedef std::vector<BcRot> EshelbianPlasticity::BcRotVec

Definition at line 563 of file EshelbianPlasticity.hpp.

◆ EntData

using EshelbianPlasticity::EntData = typedef DataForcesAndSourcesCore::EntData

Definition at line 33 of file EshelbianPlasticity.hpp.

◆ FaceUserDataOperator

Examples
NavierStokesElement.hpp.

Definition at line 36 of file EshelbianPlasticity.hpp.

◆ MatrixPtr

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

Definition at line 20 of file EshelbianPlasticity.hpp.

◆ TractionBcVec

Definition at line 574 of file EshelbianPlasticity.hpp.

◆ TractionFreeBc

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

Definition at line 565 of file EshelbianPlasticity.hpp.

◆ UserDataOperator

using EshelbianPlasticity::UserDataOperator = typedef ForcesAndSourcesCore::UserDataOperator

Definition at line 34 of file EshelbianPlasticity.hpp.

◆ VectorPtr

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

Definition at line 21 of file EshelbianPlasticity.hpp.

◆ VolUserDataOperator

Definition at line 35 of file EshelbianPlasticity.hpp.

Enumeration Type Documentation

◆ EshelbianInterfaces

Enumerator
CORE_ESHELBIAN_INTERFACE 

Definition at line 15 of file EshelbianPlasticity.hpp.

Function Documentation

◆ CGG_BubbleBase_MBTET()

MoFEMErrorCode EshelbianPlasticity::CGG_BubbleBase_MBTET ( const int  p,
const double *  N,
const double *  diffN,
FTensor::Tensor2< FTensor::PackPtr< double *, 9 >, 3, 3 > &  phi,
const int  gdim 
)

Calculate CGGT tonsorial bubble base.

See details in [14]

Parameters
ppolynomial order
Nshape functions
diffNdirevatives of shape functions
l2_basebase functions for l2 space
diff_l2_basedirevatives base functions for l2 space
phireturned base functions
gdimnumber of integration points
Returns
MoFEMErrorCode
Examples
EshelbianPlasticity.cpp.

Definition at line 20 of file CGGTonsorialBubbleBase.cpp.

23  {
25 
26  // FIXME: In this implementation symmetry of mix derivatives is not exploited
27 
28  if(p < 1)
30 
36 
37  Tensor1<double, 3> t_diff_n[4];
38  {
39  Tensor1<PackPtr<const double *, 3>, 3> t_diff_n_tmp(&diffN[0], &diffN[1],
40  &diffN[2]);
41  for (int ii = 0; ii != 4; ++ii) {
42  t_diff_n[ii](i) = t_diff_n_tmp(i);
43  ++t_diff_n_tmp;
44  }
45  }
46 
47  Tensor1<double, 3> t_diff_ksi[3];
48  for (int ii = 0; ii != 3; ++ii)
49  t_diff_ksi[ii](i) = t_diff_n[ii + 1](i) - t_diff_n[0](i);
50 
51  int lp = p >= 2 ? p - 2 + 1 : 0;
52  VectorDouble l[3] = {VectorDouble(lp + 1), VectorDouble(lp + 1),
53  VectorDouble(lp + 1)};
54  MatrixDouble diff_l[3] = {MatrixDouble(3, lp + 1), MatrixDouble(3, lp + 1),
55  MatrixDouble(3, lp + 1)};
56  MatrixDouble diff2_l[3] = {MatrixDouble(9, lp + 1), MatrixDouble(9, lp + 1),
57  MatrixDouble(9, lp + 1)};
58 
59  for (int ii = 0; ii != 3;++ii)
60  diff2_l[ii].clear();
61 
62  for (int gg = 0; gg != gdim; ++gg) {
63 
64  const int node_shift = gg * 4;
65 
66  for (int ii = 0; ii != 3; ++ii) {
67 
68  auto &t_diff_ksi_ii = t_diff_ksi[ii];
69  auto &l_ii = l[ii];
70  auto &diff_l_ii = diff_l[ii];
71  auto &diff2_l_ii = diff2_l[ii];
72 
73  double ksi_ii = N[node_shift + ii + 1] - N[node_shift + 0];
74 
75  CHKERR Legendre_polynomials(lp, ksi_ii, &t_diff_ksi_ii(0),
76  &*l_ii.data().begin(),
77  &*diff_l_ii.data().begin(), 3);
78 
79  for (int l = 1; l < lp; ++l) {
80  const double a = ((2 * (double)l + 1) / ((double)l + 1));
81  const double b = ((double)l / ((double)l + 1));
82  for (int d0 = 0; d0 != 3; ++d0)
83  for (int d1 = 0; d1 != 3; ++d1) {
84  const int r = 3 * d0 + d1;
85  diff2_l_ii(r, l + 1) = a * (t_diff_ksi_ii(d0) * diff_l_ii(d1, l) +
86  t_diff_ksi_ii(d1) * diff_l_ii(d0, l) +
87  ksi_ii * diff2_l_ii(r, l)) -
88  b * diff2_l_ii(r, l - 1);
89  }
90  }
91 
92  }
93 
94  const double n[] = {N[node_shift + 0], N[node_shift + 1], N[node_shift + 2],
95  N[node_shift + 3]};
96 
98  Tensor3<double, 3, 3, 3> t_bk_diff;
99  const int tab[4][4] = {
100  {1, 2, 3, 0}, {2, 3, 0, 1}, {3, 0, 1, 2}, {0, 1, 2, 3}};
101  t_bk(i, j) = 0;
102  t_bk_diff(i, j, k) = 0;
103  for (int ii = 0; ii != 3; ++ii) {
104  const int i0 = tab[ii][0];
105  const int i1 = tab[ii][1];
106  const int i2 = tab[ii][2];
107  const int i3 = tab[ii][3];
108  auto &t_diff_n_i0 = t_diff_n[i0];
109  auto &t_diff_n_i1 = t_diff_n[i1];
110  auto &t_diff_n_i2 = t_diff_n[i2];
111  auto &t_diff_n_i3 = t_diff_n[i3];
113  t_k(i, j) = t_diff_n_i3(i) * t_diff_n_i3(j);
114  const double b = n[i0] * n[i1] * n[i2];
115  t_bk(i, j) += b * t_k(i, j);
116  Tensor1<double, 3> t_diff_b;
117  t_diff_b(i) = t_diff_n_i0(i) * n[i1] * n[i2] +
118  t_diff_n_i1(i) * n[i0] * n[i2] +
119  t_diff_n_i2(i) * n[i0] * n[i1];
120  t_bk_diff(i, j, k) += t_k(i, j) * t_diff_b(k);
121  }
122 
123  int zz = 0;
124  for (int o = p - 2 + 1; o <= p - 2 + 1; ++o) {
125 
126  for (int ii = 0; ii <= o; ++ii)
127  for (int jj = 0; (ii + jj) <= o; ++jj) {
128 
129  const int kk = o - ii - jj;
130 
131  auto get_diff_l = [&](const int y, const int i) {
132  return Tensor1<double, 3>(diff_l[y](0, i), diff_l[y](1, i),
133  diff_l[y](2, i));
134  };
135  auto get_diff2_l = [&](const int y, const int i) {
136  return Tensor2<double, 3, 3>(
137  diff2_l[y](0, i), diff2_l[y](1, i), diff2_l[y](2, i),
138  diff2_l[y](3, i), diff2_l[y](4, i), diff2_l[y](5, i),
139  diff2_l[y](6, i), diff2_l[y](7, i), diff2_l[y](8, i));
140  };
141 
142  auto l_i = l[0][ii];
143  auto t_diff_i = get_diff_l(0, ii);
144  auto t_diff2_i = get_diff2_l(0, ii);
145  auto l_j = l[1][jj];
146  auto t_diff_j = get_diff_l(1, jj);
147  auto t_diff2_j = get_diff2_l(1, jj);
148  auto l_k = l[2][kk];
149  auto t_diff_k = get_diff_l(2, kk);
150  auto t_diff2_k = get_diff2_l(2, kk);
151 
152  Tensor1<double, 3> t_diff_l2;
153  t_diff_l2(i) = t_diff_i(i) * l_j * l_k + t_diff_j(i) * l_i * l_k +
154  t_diff_k(i) * l_i * l_j;
155  Tensor2<double, 3, 3> t_diff2_l2;
156  t_diff2_l2(i, j) =
157  t_diff2_i(i, j) * l_j * l_k + t_diff_i(i) * t_diff_j(j) * l_k +
158  t_diff_i(i) * l_j * t_diff_k(j) +
159 
160  t_diff2_j(i, j) * l_i * l_k + t_diff_j(i) * t_diff_i(j) * l_k +
161  t_diff_j(i) * l_i * t_diff_k(j) +
162 
163  t_diff2_k(i, j) * l_i * l_j + t_diff_k(i) * t_diff_i(j) * l_j +
164  t_diff_k(i) * l_i * t_diff_j(j);
165 
166  for (int dd = 0; dd != 3; ++dd) {
167 
168  Tensor2<double, 3, 3> t_axial_diff;
169  t_axial_diff(i, j) = 0;
170  for (int mm = 0; mm != 3; ++mm)
171  t_axial_diff(dd, mm) = t_diff_l2(mm);
172 
173  Tensor3<double, 3, 3, 3> t_A_diff;
174  t_A_diff(i, j, k) = levi_civita(i, j, m) * t_axial_diff(m, k);
175  Tensor2<double, 3, 3> t_curl_A;
176  t_curl_A(i, j) = levi_civita(j, m, f) * t_A_diff(i, f, m);
177  Tensor3<double, 3, 3, 3> t_curl_A_bK_diff;
178  t_curl_A_bK_diff(i, j, k) = t_curl_A(i, m) * t_bk_diff(m, j, k);
179 
180  Tensor3<double, 3, 3, 3> t_axial_diff2;
181  t_axial_diff2(i, j, k) = 0;
182  for (int mm = 0; mm != 3; ++mm)
183  for (int nn = 0; nn != 3; ++nn)
184  t_axial_diff2(dd, mm, nn) = t_diff2_l2(mm, nn);
185  Tensor4<double, 3, 3, 3, 3> t_A_diff2;
186  t_A_diff2(i, j, k, f) =
187  levi_civita(i, j, m) * t_axial_diff2(m, k, f);
188  Tensor3<double, 3, 3, 3> t_curl_A_diff2;
189  t_curl_A_diff2(i, j, k) =
190  levi_civita(j, m, f) * t_A_diff2(i, f, m, k);
191  Tensor3<double, 3, 3, 3> t_curl_A_diff2_bK;
192  t_curl_A_diff2_bK(i, j, k) = t_curl_A_diff2(i, m, k) * t_bk(m, j);
193 
194  t_phi(i, j) = levi_civita(j, m, f) * (t_curl_A_bK_diff(i, f, m) +
195  t_curl_A_diff2_bK(i, f, m));
196 
197  ++t_phi;
198  ++zz;
199  }
200  }
201  }
202  if(zz != NBVOLUMETET_CCG_BUBBLE(p))
203  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
204  "Wrong number of base functions %d != %d", zz,
206 
207  }
208 
210 }
static Index< 'm', 3 > m
static Index< 'n', 3 > n
static Index< 'o', 3 > o
static Index< 'p', 3 > p
#define NBVOLUMETET_CCG_BUBBLE(P)
Bubble function for CGG H div space.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:123
#define CHKERR
Inline error check.
Definition: definitions.h:604
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
PetscErrorCode Legendre_polynomials(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
Calculate Legendre approximation basis.
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'i', 3 > i
FTensor::Index< 'k', 3 > k
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
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
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:74
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:76
const double r
rate factor
const int N
Definition: speed_test.cpp:3

◆ d_f()

double EshelbianPlasticity::d_f ( const double  v)
Examples
EshelbianOperators.cpp.

Definition at line 23 of file EshelbianOperators.cpp.

23 { return exp(v); }

◆ dd_f()

double EshelbianPlasticity::dd_f ( const double  v)
Examples
EshelbianOperators.cpp.

Definition at line 24 of file EshelbianOperators.cpp.

24 { return exp(v); }

◆ f()

double EshelbianPlasticity::f ( const double  v)
Examples
EshelbianOperators.cpp.

Definition at line 22 of file EshelbianOperators.cpp.

22 { return exp(v); }

◆ get_diff_rotation_form_vector()

template<typename T >
static FTensor::Tensor3<typename TensorTypeExtractor<T>::Type, 3, 3, 3> EshelbianPlasticity::get_diff_rotation_form_vector ( FTensor::Tensor1< T, 3 > &  t_omega)
static
Examples
EshelbianOperators.cpp.

Definition at line 66 of file EshelbianOperators.cpp.

67  {
68 
73 
74  const typename TensorTypeExtractor<T>::Type angle =
75  sqrt(t_omega(i) * t_omega(i));
76 
77  constexpr double tol = 1e-18;
78  if (std::abs(angle) < tol) {
80  auto t_R = get_rotation_form_vector(t_omega);
81  t_diff_R(i, j, k) = FTensor::levi_civita<double>(i, l, k) * t_R(l, j);
82  return t_diff_R;
83  }
84 
87  t_Omega(i, j) = FTensor::levi_civita<double>(i, j, k) * t_omega(k);
88 
89  const typename TensorTypeExtractor<T>::Type angle2 = angle * angle;
90  const typename TensorTypeExtractor<T>::Type ss = sin(angle);
91  const typename TensorTypeExtractor<T>::Type a = ss / angle;
92  const typename TensorTypeExtractor<T>::Type ss_2 = sin(angle / 2.);
93  const typename TensorTypeExtractor<T>::Type b = 2. * ss_2 * ss_2 / angle2;
94  t_diff_R(i, j, k) = 0;
95  t_diff_R(i, j, k) = a * FTensor::levi_civita<double>(i, j, k);
96  t_diff_R(i, j, k) +=
97  b * (t_Omega(i, l) * FTensor::levi_civita<double>(l, j, k) +
98  FTensor::levi_civita<double>(i, l, k) * t_Omega(l, j));
99 
100  const typename TensorTypeExtractor<T>::Type cc = cos(angle);
101  const typename TensorTypeExtractor<T>::Type cc_2 = cos(angle / 2.);
102  const typename TensorTypeExtractor<T>::Type diff_a =
103  (angle * cc - ss) / angle2;
104 
105  const typename TensorTypeExtractor<T>::Type diff_b =
106  (2. * angle * ss_2 * cc_2 - 4. * ss_2 * ss_2) / (angle2 * angle);
107  t_diff_R(i, j, k) += diff_a * t_Omega(i, j) * (t_omega(k) / angle);
108  t_diff_R(i, j, k) +=
109  diff_b * t_Omega(i, l) * t_Omega(l, j) * (t_omega(k) / angle);
110 
111  return t_diff_R;
112 }
double tol
static FTensor::Tensor2< typename TensorTypeExtractor< T >::Type, 3, 3 > get_rotation_form_vector(FTensor::Tensor1< T, 3 > &t_omega)

◆ get_eigen_val_and_proj_lapack()

template<typename T1 , typename T2 , typename T3 >
MoFEMErrorCode EshelbianPlasticity::get_eigen_val_and_proj_lapack ( T1 &  X,
T2 &  eig,
T3 &  eig_vec 
)

Definition at line 115 of file EshelbianOperators.cpp.

116  {
118  for (int ii = 0; ii != 3; ii++)
119  for (int jj = 0; jj != 3; jj++)
120  eig_vec(ii, jj) = X(ii, jj);
121 
122  int n = 3;
123  int lda = 3;
124  int lwork = (3 + 2) * 3;
125  std::array<double, (3 + 2) * 3> work;
126 
127  if (lapack_dsyev('V', 'U', n, &(eig_vec(0, 0)), lda, &eig(0), work.data(),
128  lwork) > 0)
129  SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA,
130  "The algorithm failed to compute eigenvalues.");
132 }
@ MOFEM_INVALID_DATA
Definition: definitions.h:128
static __CLPK_integer lapack_dsyev(char jobz, char uplo, __CLPK_integer n, __CLPK_doublereal *a, __CLPK_integer lda, __CLPK_doublereal *w, __CLPK_doublereal *work, __CLPK_integer lwork)
Definition: lapack_wrap.h:273

◆ get_rotation_form_vector()

template<typename T >
static FTensor::Tensor2<typename TensorTypeExtractor<T>::Type, 3, 3> EshelbianPlasticity::get_rotation_form_vector ( FTensor::Tensor1< T, 3 > &  t_omega)
static
Examples
EshelbianOperators.cpp.

Definition at line 35 of file EshelbianOperators.cpp.

36  {
41  constexpr auto t_kd = FTensor::Kronecker_Delta<int>();
42  t_R(i, j) = t_kd(i, j);
43 
44  const typename TensorTypeExtractor<T>::Type angle =
45  sqrt(t_omega(i) * t_omega(i));
46 
47  constexpr double tol = 1e-18;
48  if (std::abs(angle) < tol) {
49  return t_R;
50  }
51 
53  t_Omega(i, j) = FTensor::levi_civita<double>(i, j, k) * t_omega(k);
54  const typename TensorTypeExtractor<T>::Type a = sin(angle) / angle;
55  const typename TensorTypeExtractor<T>::Type ss_2 = sin(angle / 2.);
56  const typename TensorTypeExtractor<T>::Type b =
57  2. * ss_2 * ss_2 / (angle * angle);
58  t_R(i, j) += a * t_Omega(i, j);
59  t_R(i, j) += b * t_Omega(i, k) * t_Omega(k, j);
60 
61  return t_R;
62 }
Kronecker Delta class.

◆ get_uniq_nb()

template<typename T >
int EshelbianPlasticity::get_uniq_nb ( T ec)
Examples
EshelbianOperators.cpp.

Definition at line 139 of file EshelbianOperators.cpp.

139  {
140  return distance(&ec(0), unique(&ec(0), &ec(0) + 3, is_eq));
141 }
bool is_eq(const double &a, const double &b)

◆ getFTensor3FromMat()

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

Definition at line 24 of file EshelbianPlasticity.hpp.

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

◆ is_eq()

bool EshelbianPlasticity::is_eq ( const double &  a,
const double &  b 
)

Definition at line 134 of file EshelbianOperators.cpp.

134  {
135  constexpr double eps = 1e-8;
136  return abs(a - b) < eps;
137 }

◆ sort_eigen_vals()

template<typename T1 , typename T2 >
void EshelbianPlasticity::sort_eigen_vals ( T1 &  eig,
T2 &  eigen_vec 
)
Examples
EshelbianOperators.cpp.

Definition at line 144 of file EshelbianOperators.cpp.

144  {
147  if (is_eq(eig(0), eig(1))) {
148  FTensor::Tensor2<double, 3, 3> eigen_vec_c{
149  eigen_vec(0, 0), eigen_vec(0, 1), eigen_vec(0, 2),
150  eigen_vec(2, 0), eigen_vec(2, 1), eigen_vec(2, 2),
151  eigen_vec(1, 0), eigen_vec(1, 1), eigen_vec(1, 2)};
152  FTensor::Tensor1<double, 3> eig_c{eig(0), eig(2), eig(1)};
153  eig(i) = eig_c(i);
154  eigen_vec(i, j) = eigen_vec_c(i, j);
155  } else {
156  FTensor::Tensor2<double, 3, 3> eigen_vec_c{
157  eigen_vec(1, 0), eigen_vec(1, 1), eigen_vec(1, 2),
158  eigen_vec(0, 0), eigen_vec(0, 1), eigen_vec(0, 2),
159  eigen_vec(2, 0), eigen_vec(2, 1), eigen_vec(2, 2)};
160  FTensor::Tensor1<double, 3> eig_c{eig(1), eig(0), eig(2)};
161  eig(i) = eig_c(i);
162  eigen_vec(i, j) = eigen_vec_c(i, j);
163  }
164 }

Variable Documentation

◆ IDD_MOFEMEshelbianCrackInterface

const MOFEMuuid EshelbianPlasticity::IDD_MOFEMEshelbianCrackInterface
static
Initial value:
=
std::bitset< BITINTERFACEUID_SIZE > BitIntefaceId
Definition: Types.hpp:56
Examples
EshelbianPlasticity.cpp.

Definition at line 17 of file EshelbianPlasticity.hpp.