v0.14.0
Loading...
Searching...
No Matches
Classes | Typedefs | Enumerations | Functions | Variables
EshelbianPlasticity Namespace Reference

Classes

struct  BcDisp
 
struct  BcRot
 
struct  CGGUserPolynomialBase
 
struct  ContactTree
 
struct  DataAtIntegrationPts
 
struct  EshelbianCore
 
struct  FaceRule
 
struct  FeTractionBc
 
struct  HMHHencky
 
struct  HMHPMooneyRivlinWriggersEq63
 
struct  HMHStVenantKirchhoff
 
struct  isEq
 
struct  OpAssembleBasic
 
struct  OpAssembleFace
 
struct  OpAssembleVolume
 
struct  OpCalculateRotationAndSpatialGradient
 
struct  OpCalculateStrainEnergy
 
struct  OpConstrainBoundaryLhs_dP
 
struct  OpConstrainBoundaryLhs_dU
 
struct  OpConstrainBoundaryRhs
 
struct  OpDispBc
 
struct  OpDispBc_dx
 
struct  OpHMHH
 
struct  OpJacobian
 
struct  OpLoopSideGetDataForSideEle
 
struct  OpMoveNode
 
struct  OpPostProcDataStructure
 
struct  OpRotationBc
 
struct  OpRotationBc_dx
 
struct  OpSpatialConsistency_dBubble_domega
 
struct  OpSpatialConsistency_dP_domega
 
struct  OpSpatialConsistencyBubble
 
struct  OpSpatialConsistencyDivTerm
 
struct  OpSpatialConsistencyP
 
struct  OpSpatialEquilibrium
 
struct  OpSpatialEquilibrium_dw_dP
 
struct  OpSpatialEquilibrium_dw_dw
 
struct  OpSpatialPhysical
 
struct  OpSpatialPhysical_du_dBubble
 
struct  OpSpatialPhysical_du_domega
 
struct  OpSpatialPhysical_du_dP
 
struct  OpSpatialPhysical_du_du
 
struct  OpSpatialPhysical_du_dx
 
struct  OpSpatialPrj
 
struct  OpSpatialPrj_dx_dw
 
struct  OpSpatialPrj_dx_dx
 
struct  OpSpatialRotation
 
struct  OpSpatialRotation_domega_dBubble
 
struct  OpSpatialRotation_domega_domega
 
struct  OpSpatialRotation_domega_dP
 
struct  OpTractionBc
 
struct  OpTreeSearch
 
struct  PhysicalEquations
 
struct  TensorTypeExtractor
 
struct  TensorTypeExtractor< FTensor::PackPtr< T, I > >
 
struct  TractionBc
 
struct  VolRule
 Set integration rule on element. More...
 
struct  VolUserDataOperatorStabAssembly
 

Typedefs

using MatrixPtr = boost::shared_ptr< MatrixDouble >
 
using VectorPtr = boost::shared_ptr< VectorDouble >
 
using EntData = EntitiesFieldData::EntData
 
using UserDataOperator = ForcesAndSourcesCore::UserDataOperator
 
using VolUserDataOperator = VolumeElementForcesAndSourcesCore::UserDataOperator
 
using FaceUserDataOperator = FaceElementForcesAndSourcesCore::UserDataOperator
 
typedef std::vector< BcDispBcDispVec
 
typedef std::vector< BcRotBcRotVec
 
typedef std::vector< RangeTractionFreeBc
 
typedef std::vector< TractionBcTractionBcVec
 

Enumerations

enum  RotSelector { SMALL_ROT , MODERATE_ROT , LARGE_ROT }
 
enum  StretchSelector { LINEAR , LOG }
 

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...
 
auto diff_deviator (FTensor::Ddg< double, 3, 3 > &&t_diff_stress)
 
auto diff_tensor ()
 
auto symm_L_tensor ()
 
template<typename T >
auto get_rotation_form_vector (FTensor::Tensor1< T, 3 > &t_omega, RotSelector rotSelector=LARGE_ROT)
 
template<typename T >
auto get_diff_rotation_form_vector (FTensor::Tensor1< T, 3 > &t_omega, RotSelector rotSelector=LARGE_ROT)
 
template<typename T >
auto get_diff2_rotation_form_vector (FTensor::Tensor1< T, 3 > &t_omega, RotSelector rotSelector=LARGE_ROT)
 
auto is_eq (const double &a, const double &b)
 
template<int DIM>
auto get_uniq_nb (double *ptr)
 
template<int DIM>
auto sort_eigen_vals (FTensor::Tensor1< double, DIM > &eig, FTensor::Tensor2< double, DIM, DIM > &eigen_vec)
 

Variables

static constexpr auto size_symm = (3 * (3 + 1)) / 2
 
enum RotSelector EshelbianCore
 

Typedef Documentation

◆ BcDispVec

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

Definition at line 298 of file EshelbianPlasticity.hpp.

◆ BcRotVec

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

Definition at line 307 of file EshelbianPlasticity.hpp.

◆ EntData

using EshelbianPlasticity::EntData = typedef EntitiesFieldData::EntData

Definition at line 26 of file EshelbianPlasticity.hpp.

◆ FaceUserDataOperator

using EshelbianPlasticity::FaceUserDataOperator = typedef FaceElementForcesAndSourcesCore::UserDataOperator

Definition at line 29 of file EshelbianPlasticity.hpp.

◆ MatrixPtr

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

Definition at line 23 of file EshelbianPlasticity.hpp.

◆ TractionBcVec

Definition at line 318 of file EshelbianPlasticity.hpp.

◆ TractionFreeBc

Definition at line 309 of file EshelbianPlasticity.hpp.

◆ UserDataOperator

using EshelbianPlasticity::UserDataOperator = typedef ForcesAndSourcesCore::UserDataOperator

Definition at line 27 of file EshelbianPlasticity.hpp.

◆ VectorPtr

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

Definition at line 24 of file EshelbianPlasticity.hpp.

◆ VolUserDataOperator

using EshelbianPlasticity::VolUserDataOperator = typedef VolumeElementForcesAndSourcesCore::UserDataOperator

Definition at line 28 of file EshelbianPlasticity.hpp.

Enumeration Type Documentation

◆ RotSelector

Enumerator
SMALL_ROT 
MODERATE_ROT 
LARGE_ROT 
Examples
EshelbianOperators.cpp, and EshelbianPlasticity.cpp.

Definition at line 20 of file EshelbianPlasticity.hpp.

◆ StretchSelector

Enumerator
LINEAR 
LOG 
Examples
EshelbianPlasticity.cpp.

Definition at line 21 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 [20]

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

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 const double n[] = {N[node_shift + 0], N[node_shift + 1], N[node_shift + 2],
94 N[node_shift + 3]};
95
98 const int tab[4][4] = {
99 {1, 2, 3, 0}, {2, 3, 0, 1}, {3, 0, 1, 2}, {0, 1, 2, 3}};
100 t_bk(i, j) = 0;
101 t_bk_diff(i, j, k) = 0;
102 for (int ii = 0; ii != 3; ++ii) {
103 const int i0 = tab[ii][0];
104 const int i1 = tab[ii][1];
105 const int i2 = tab[ii][2];
106 const int i3 = tab[ii][3];
107 auto &t_diff_n_i0 = t_diff_n[i0];
108 auto &t_diff_n_i1 = t_diff_n[i1];
109 auto &t_diff_n_i2 = t_diff_n[i2];
110 auto &t_diff_n_i3 = t_diff_n[i3];
112 t_k(i, j) = t_diff_n_i3(i) * t_diff_n_i3(j);
113 const double b = n[i0] * n[i1] * n[i2];
114 t_bk(i, j) += b * t_k(i, j);
115 Tensor1<double, 3> t_diff_b;
116 t_diff_b(i) = t_diff_n_i0(i) * n[i1] * n[i2] +
117 t_diff_n_i1(i) * n[i0] * n[i2] +
118 t_diff_n_i2(i) * n[i0] * n[i1];
119 t_bk_diff(i, j, k) += t_k(i, j) * t_diff_b(k);
120 }
121
122 int zz = 0;
123 for (int o = p - 2 + 1; o <= p - 2 + 1; ++o) {
124
125 for (int ii = 0; ii <= o; ++ii)
126 for (int jj = 0; (ii + jj) <= o; ++jj) {
127
128 const int kk = o - ii - jj;
129
130 auto get_diff_l = [&](const int y, const int i) {
131 return Tensor1<double, 3>(diff_l[y](0, i), diff_l[y](1, i),
132 diff_l[y](2, i));
133 };
134 auto get_diff2_l = [&](const int y, const int i) {
136 diff2_l[y](0, i), diff2_l[y](1, i), diff2_l[y](2, i),
137 diff2_l[y](3, i), diff2_l[y](4, i), diff2_l[y](5, i),
138 diff2_l[y](6, i), diff2_l[y](7, i), diff2_l[y](8, i));
139 };
140
141 auto l_i = l[0][ii];
142 auto t_diff_i = get_diff_l(0, ii);
143 auto t_diff2_i = get_diff2_l(0, ii);
144 auto l_j = l[1][jj];
145 auto t_diff_j = get_diff_l(1, jj);
146 auto t_diff2_j = get_diff2_l(1, jj);
147 auto l_k = l[2][kk];
148 auto t_diff_k = get_diff_l(2, kk);
149 auto t_diff2_k = get_diff2_l(2, kk);
150
151 Tensor1<double, 3> t_diff_l2;
152 t_diff_l2(i) = t_diff_i(i) * l_j * l_k + t_diff_j(i) * l_i * l_k +
153 t_diff_k(i) * l_i * l_j;
154 Tensor2<double, 3, 3> t_diff2_l2;
155 t_diff2_l2(i, j) =
156 t_diff2_i(i, j) * l_j * l_k + t_diff_i(i) * t_diff_j(j) * l_k +
157 t_diff_i(i) * l_j * t_diff_k(j) +
158
159 t_diff2_j(i, j) * l_i * l_k + t_diff_j(i) * t_diff_i(j) * l_k +
160 t_diff_j(i) * l_i * t_diff_k(j) +
161
162 t_diff2_k(i, j) * l_i * l_j + t_diff_k(i) * t_diff_i(j) * l_j +
163 t_diff_k(i) * l_i * t_diff_j(j);
164
165 for (int dd = 0; dd != 3; ++dd) {
166
167 Tensor2<double, 3, 3> t_axial_diff;
168 t_axial_diff(i, j) = 0;
169 for (int mm = 0; mm != 3; ++mm)
170 t_axial_diff(dd, mm) = t_diff_l2(mm);
171
173 t_A_diff(i, j, k) = levi_civita(i, j, m) * t_axial_diff(m, k);
174 Tensor2<double, 3, 3> t_curl_A;
175 t_curl_A(i, j) = levi_civita(j, m, f) * t_A_diff(i, f, m);
176 Tensor3<double, 3, 3, 3> t_curl_A_bK_diff;
177 t_curl_A_bK_diff(i, j, k) = t_curl_A(i, m) * t_bk_diff(m, j, k);
178
179 Tensor3<double, 3, 3, 3> t_axial_diff2;
180 t_axial_diff2(i, j, k) = 0;
181 for (int mm = 0; mm != 3; ++mm)
182 for (int nn = 0; nn != 3; ++nn)
183 t_axial_diff2(dd, mm, nn) = t_diff2_l2(mm, nn);
185 t_A_diff2(i, j, k, f) =
186 levi_civita(i, j, m) * t_axial_diff2(m, k, f);
187 Tensor3<double, 3, 3, 3> t_curl_A_diff2;
188 t_curl_A_diff2(i, j, k) =
189 levi_civita(j, m, f) * t_A_diff2(i, f, m, k);
190 Tensor3<double, 3, 3, 3> t_curl_A_diff2_bK;
191 t_curl_A_diff2_bK(i, j, k) = t_curl_A_diff2(i, m, k) * t_bk(m, j);
192
193 t_phi(i, j) = levi_civita(j, m, f) * (t_curl_A_bK_diff(i, f, m) +
194 t_curl_A_diff2_bK(i, f, m));
195
196 ++t_phi;
197 ++zz;
198 }
199 }
200 }
201 if (zz != NBVOLUMETET_CCG_BUBBLE(p))
202 SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
203 "Wrong number of base functions %d != %d", zz,
205 }
206
208}
static Index< 'o', 3 > o
static Index< 'p', 3 > p
#define NBVOLUMETET_CCG_BUBBLE(P)
Bubble function for CGG H div space.
constexpr double a
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
#define CHKERR
Inline error check.
Definition: definitions.h:535
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
FTensor::Index< 'n', SPACE_DIM > n
FTensor::Index< 'm', SPACE_DIM > m
PetscErrorCode Legendre_polynomials(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
Calculate Legendre approximation basis.
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
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
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
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:77
UBlasVector< double > VectorDouble
Definition: Types.hpp:68
int r
Definition: sdf.py:8
const int N
Definition: speed_test.cpp:3

◆ diff_deviator()

auto EshelbianPlasticity::diff_deviator ( FTensor::Ddg< double, 3, 3 > &&  t_diff_stress)
inline
Examples
EshelbianOperators.cpp.

Definition at line 20 of file EshelbianOperators.cpp.

20 {
25
26 FTensor::Ddg<double, 3, 3> t_diff_deviator;
27 t_diff_deviator(i, j, k, l) = t_diff_stress(i, j, k, l);
28
29 constexpr double third = boost::math::constants::third<double>();
30
31 t_diff_deviator(0, 0, 0, 0) -= third;
32 t_diff_deviator(0, 0, 1, 1) -= third;
33
34 t_diff_deviator(1, 1, 0, 0) -= third;
35 t_diff_deviator(1, 1, 1, 1) -= third;
36
37 t_diff_deviator(2, 2, 0, 0) -= third;
38 t_diff_deviator(2, 2, 1, 1) -= third;
39
40 t_diff_deviator(0, 0, 2, 2) -= third;
41 t_diff_deviator(1, 1, 2, 2) -= third;
42 t_diff_deviator(2, 2, 2, 2) -= third;
43
44 return t_diff_deviator;
45}
constexpr double third

◆ diff_tensor()

auto EshelbianPlasticity::diff_tensor ( )
inline
Examples
EshelbianOperators.cpp.

Definition at line 49 of file EshelbianOperators.cpp.

49 {
56 t_diff(i, j, k, l) = (t_kd(i, k) ^ t_kd(j, l)) / 4.;
57 return t_diff;
58};
Kronecker Delta class symmetric.
constexpr auto t_kd

◆ get_diff2_rotation_form_vector()

template<typename T >
auto EshelbianPlasticity::get_diff2_rotation_form_vector ( FTensor::Tensor1< T, 3 > &  t_omega,
RotSelector  rotSelector = LARGE_ROT 
)
Examples
EshelbianOperators.cpp.

Definition at line 167 of file EshelbianOperators.cpp.

168 {
169
170 using D = typename TensorTypeExtractor<T>::Type;
171
174
175 constexpr double eps = 1e-10;
176 for (int l = 0; l != 3; ++l) {
178 t_omega_c(i) = t_omega(i);
179 t_omega_c(l) += std::complex<double>(0, eps);
180 auto t_diff_R_c = get_diff_rotation_form_vector(t_omega_c, rotSelector);
181 for (int i = 0; i != 3; ++i) {
182 for (int j = 0; j != 3; ++j) {
183 for (int k = 0; k != 3; ++k) {
184 t_diff2_R(i, j, k, l) = t_diff_R_c(i, j, k).imag() / eps;
185 }
186 }
187 }
188 }
189
190 return t_diff2_R;
191}
static const double eps
double D
auto get_diff_rotation_form_vector(FTensor::Tensor1< T, 3 > &t_omega, RotSelector rotSelector=LARGE_ROT)

◆ get_diff_rotation_form_vector()

template<typename T >
auto EshelbianPlasticity::get_diff_rotation_form_vector ( FTensor::Tensor1< T, 3 > &  t_omega,
RotSelector  rotSelector = LARGE_ROT 
)
Examples
EshelbianOperators.cpp.

Definition at line 119 of file EshelbianOperators.cpp.

120 {
121
122 using D = typename TensorTypeExtractor<T>::Type;
123
128
130
131 const auto angle =
132 sqrt(t_omega(i) * t_omega(i)) + std::numeric_limits<double>::epsilon();
133 if (rotSelector == SMALL_ROT) {
134 t_diff_R(i, j, k) = FTensor::levi_civita<double>(i, j, k);
135 return t_diff_R;
136 }
137
138 const auto ss = sin(angle);
139 const auto a = ss / angle;
140 t_diff_R(i, j, k) = a * FTensor::levi_civita<double>(i, j, k);
141
143 t_Omega(i, j) = FTensor::levi_civita<double>(i, j, k) * t_omega(k);
144
145 const auto angle2 = angle * angle;
146 const auto cc = cos(angle);
147 const auto diff_a = (angle * cc - ss) / angle2;
148 t_diff_R(i, j, k) += diff_a * t_Omega(i, j) * (t_omega(k) / angle);
149 if (rotSelector == MODERATE_ROT)
150 return t_diff_R;
151
152 const auto ss_2 = sin(angle / 2.);
153 const auto cc_2 = cos(angle / 2.);
154 const auto b = 2. * ss_2 * ss_2 / angle2;
155 t_diff_R(i, j, k) +=
156 b * (t_Omega(i, l) * FTensor::levi_civita<double>(l, j, k) +
157 FTensor::levi_civita<double>(i, l, k) * t_Omega(l, j));
158 const auto diff_b =
159 (2. * angle * ss_2 * cc_2 - 4. * ss_2 * ss_2) / (angle2 * angle);
160 t_diff_R(i, j, k) +=
161 diff_b * t_Omega(i, l) * t_Omega(l, j) * (t_omega(k) / angle);
162
163 return t_diff_R;
164}

◆ get_rotation_form_vector()

template<typename T >
auto EshelbianPlasticity::get_rotation_form_vector ( FTensor::Tensor1< T, 3 > &  t_omega,
RotSelector  rotSelector = LARGE_ROT 
)
Examples
EshelbianOperators.cpp.

Definition at line 83 of file EshelbianOperators.cpp.

84 {
85
86 using D = typename TensorTypeExtractor<T>::Type;
87
92
93 constexpr auto t_kd = FTensor::Kronecker_Delta<int>();
94 t_R(i, j) = t_kd(i, j);
95
97 t_Omega(i, j) = FTensor::levi_civita<double>(i, j, k) * t_omega(k);
98
99 const auto angle =
100 sqrt(t_omega(i) * t_omega(i)) + std::numeric_limits<double>::epsilon();
101 if (rotSelector == SMALL_ROT) {
102 t_R(i, j) += t_Omega(i, j);
103 return t_R;
104 }
105
106 const auto a = sin(angle) / angle;
107 t_R(i, j) += a * t_Omega(i, j);
108 if (rotSelector == MODERATE_ROT)
109 return t_R;
110
111 const auto ss_2 = sin(angle / 2.);
112 const auto b = 2. * ss_2 * ss_2 / (angle * angle);
113 t_R(i, j) += b * t_Omega(i, k) * t_Omega(k, j);
114
115 return t_R;
116}
Kronecker Delta class.

◆ get_uniq_nb()

template<int DIM>
auto EshelbianPlasticity::get_uniq_nb ( double ptr)
inline
Examples
EshelbianOperators.cpp.

Definition at line 207 of file EshelbianOperators.cpp.

207 {
208 std::array<double, DIM> tmp;
209 std::copy(ptr, &ptr[DIM], tmp.begin());
210 std::sort(tmp.begin(), tmp.end());
211 isEq::absMax = std::max(std::abs(tmp[0]), std::abs(tmp[DIM - 1]));
212 constexpr double eps = std::numeric_limits<float>::epsilon();
213 isEq::absMax = std::max(isEq::absMax, static_cast<double>(eps));
214 return std::distance(tmp.begin(), std::unique(tmp.begin(), tmp.end(), is_eq));
215}

◆ is_eq()

auto EshelbianPlasticity::is_eq ( const double a,
const double b 
)
inline
Examples
EshelbianOperators.cpp.

Definition at line 203 of file EshelbianOperators.cpp.

203 {
204 return isEq::check(a, b);
205};

◆ sort_eigen_vals()

template<int DIM>
auto EshelbianPlasticity::sort_eigen_vals ( FTensor::Tensor1< double, DIM > &  eig,
FTensor::Tensor2< double, DIM, DIM > &  eigen_vec 
)
inline
Examples
EshelbianOperators.cpp.

Definition at line 218 of file EshelbianOperators.cpp.

219 {
220
221 isEq::absMax =
222 std::max(std::max(std::abs(eig(0)), std::abs(eig(1))), std::abs(eig(2)));
223
224 int i = 0, j = 1, k = 2;
225
226 if (is_eq(eig(0), eig(1))) {
227 i = 0;
228 j = 2;
229 k = 1;
230 } else if (is_eq(eig(0), eig(2))) {
231 i = 0;
232 j = 1;
233 k = 2;
234 } else if (is_eq(eig(1), eig(2))) {
235 i = 1;
236 j = 0;
237 k = 2;
238 }
239
241 eigen_vec(i, 0), eigen_vec(i, 1), eigen_vec(i, 2),
242
243 eigen_vec(j, 0), eigen_vec(j, 1), eigen_vec(j, 2),
244
245 eigen_vec(k, 0), eigen_vec(k, 1), eigen_vec(k, 2)};
246
247 FTensor::Tensor1<double, 3> eig_c{eig(i), eig(j), eig(k)};
248
249 {
252 eig(i) = eig_c(i);
253 eigen_vec(i, j) = eigen_vec_c(i, j);
254 }
255}
auto is_eq(const double &a, const double &b)

◆ symm_L_tensor()

auto EshelbianPlasticity::symm_L_tensor ( )
inline
Examples
EshelbianOperators.cpp.

Definition at line 60 of file EshelbianOperators.cpp.

60 {
65 t_L(i, j, L) = 0;
66 t_L(0, 0, 0) = 1;
67 t_L(1, 1, 3) = 1;
68 t_L(2, 2, 5) = 1;
69 t_L(1, 0, 1) = 1;
70 t_L(2, 0, 2) = 1;
71 t_L(2, 1, 4) = 1;
72 return t_L;
73}

Variable Documentation

◆ EshelbianCore

Definition at line 55 of file EshelbianPlasticity.cpp.

◆ size_symm

constexpr auto EshelbianPlasticity::size_symm = (3 * (3 + 1)) / 2
staticconstexpr
Examples
EshelbianOperators.cpp.

Definition at line 47 of file EshelbianOperators.cpp.