v0.13.1
HenckyOps.hpp
Go to the documentation of this file.
1namespace HenckyOps {
2
3constexpr double eps = std::numeric_limits<float>::epsilon();
4
5auto f = [](double v) { return 0.5 * std::log(static_cast<long double>(v)); };
6auto d_f = [](double v) { return 0.5 / static_cast<long double>(v); };
7auto dd_f = [](double v) {
8 return -0.5 / (static_cast<long double>(v) * static_cast<long double>(v));
9};
10
11inline auto is_eq(const double &a, const double &b) {
12 return std::abs(a - b) < eps;
13};
14
15template <int DIM> inline auto get_uniq_nb(double *ptr) {
16 std::array<double, DIM> tmp;
17 std::copy(ptr, &ptr[DIM], tmp.begin());
18 std::sort(tmp.begin(), tmp.end());
19 return std::distance(tmp.begin(), std::unique(tmp.begin(), tmp.end(), is_eq));
20};
21
22template <int DIM>
25 std::array<std::pair<double, size_t>, DIM> tmp;
26 auto is_eq_pair = [](const auto &a, const auto &b) {
27 return a.first < b.first;
28 };
29
30 for (size_t i = 0; i != DIM; ++i)
31 tmp[i] = std::make_pair(eig(i), i);
32 std::sort(tmp.begin(), tmp.end(), is_eq_pair);
33
34 int i, j, k;
35 if (is_eq_pair(tmp[0], tmp[1])) {
36 i = tmp[0].second;
37 j = tmp[2].second;
38 k = tmp[1].second;
39 } else {
40 i = tmp[0].second;
41 j = tmp[1].second;
42 k = tmp[2].second;
43 }
44
46 eigen_vec(i, 0), eigen_vec(i, 1), eigen_vec(i, 2),
47
48 eigen_vec(j, 0), eigen_vec(j, 1), eigen_vec(j, 2),
49
50 eigen_vec(k, 0), eigen_vec(k, 1), eigen_vec(k, 2)};
51
52 FTensor::Tensor1<double, 3> eig_c{eig(i), eig(j), eig(k)};
53
54 {
57 eig(i) = eig_c(i);
58 eigen_vec(i, j) = eigen_vec_c(i, j);
59 }
60};
61
62struct CommonData : public boost::enable_shared_from_this<CommonData> {
63 boost::shared_ptr<MatrixDouble> matGradPtr;
64 boost::shared_ptr<MatrixDouble> matDPtr;
65 boost::shared_ptr<MatrixDouble> matLogCPlastic;
66
75
76 inline auto getMatFirstPiolaStress() {
77 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
79 }
80
81 inline auto getMatHenckyStress() {
82 return boost::shared_ptr<MatrixDouble>(shared_from_this(),
84 }
85
86 inline auto getMatLogC() {
87 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &matLogC);
88 }
89
90 inline auto getMatTangent() {
91 return boost::shared_ptr<MatrixDouble>(shared_from_this(), &matTangent);
92 }
93};
94
95template <int DIM> struct OpCalculateEigenVals : public DomainEleOp {
96
98 boost::shared_ptr<CommonData> common_data)
100 commonDataPtr(common_data) {
101 std::fill(&doEntities[MBEDGE], &doEntities[MBMAXTYPE], false);
102 }
103
106
110
111 constexpr auto t_kd = FTensor::Kronecker_Delta<int>();
112 // const size_t nb_gauss_pts = matGradPtr->size2();
113 const size_t nb_gauss_pts = getGaussPts().size2();
114
115 auto t_grad = getFTensor2FromMat<DIM, DIM>(*(commonDataPtr->matGradPtr));
116
117 commonDataPtr->matEigVal.resize(DIM, nb_gauss_pts, false);
118 commonDataPtr->matEigVec.resize(DIM * DIM, nb_gauss_pts, false);
119 auto t_eig_val = getFTensor1FromMat<DIM>(commonDataPtr->matEigVal);
120 auto t_eig_vec = getFTensor2FromMat<DIM, DIM>(commonDataPtr->matEigVec);
121
122 for (size_t gg = 0; gg != nb_gauss_pts; ++gg) {
123
128
129 F(i, j) = t_grad(i, j) + t_kd(i, j);
130 C(i, j) = F(k, i) ^ F(k, j);
131
132 for (int ii = 0; ii != DIM; ii++)
133 for (int jj = 0; jj != DIM; jj++)
134 eigen_vec(ii, jj) = C(ii, jj);
135
136 CHKERR computeEigenValuesSymmetric<DIM>(eigen_vec, eig);
137
138 // rare case when two eigen values are equal
139 auto nb_uniq = get_uniq_nb<DIM>(&eig(0));
140 if (DIM == 3 && nb_uniq == 2)
141 sort_eigen_vals<DIM>(eig, eigen_vec);
142
143 t_eig_val(i) = eig(i);
144 t_eig_vec(i, j) = eigen_vec(i, j);
145
146 ++t_grad;
147 ++t_eig_val;
148 ++t_eig_vec;
149 }
150
152 }
153
154private:
155 boost::shared_ptr<CommonData> commonDataPtr;
156};
157
158template <int DIM> struct OpCalculateLogC : public DomainEleOp {
159
160 OpCalculateLogC(const std::string field_name,
161 boost::shared_ptr<CommonData> common_data)
163 commonDataPtr(common_data) {
164 std::fill(&doEntities[MBEDGE], &doEntities[MBMAXTYPE], false);
165 }
166
169
172
173 constexpr auto t_kd = FTensor::Kronecker_Delta<int>();
174 // const size_t nb_gauss_pts = matGradPtr->size2();
175 const size_t nb_gauss_pts = getGaussPts().size2();
176 constexpr auto size_symm = (DIM * (DIM + 1)) / 2;
177 commonDataPtr->matLogC.resize(size_symm, nb_gauss_pts, false);
178
179 auto t_eig_val = getFTensor1FromMat<DIM>(commonDataPtr->matEigVal);
180 auto t_eig_vec = getFTensor2FromMat<DIM, DIM>(commonDataPtr->matEigVec);
181
182 auto t_logC = getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matLogC);
183
184 for (size_t gg = 0; gg != nb_gauss_pts; ++gg) {
185
186 // rare case when two eigen values are equal
187 auto nb_uniq = get_uniq_nb<DIM>(&(t_eig_val(0)));
188
191 eig(i) = t_eig_val(i);
192 eigen_vec(i, j) = t_eig_vec(i, j);
193 auto logC = EigenMatrix::getMat(eig, eigen_vec, f);
194 t_logC(i, j) = logC(i, j);
195
196 ++t_eig_val;
197 ++t_eig_vec;
198 ++t_logC;
199 }
200
202 }
203
204private:
205 boost::shared_ptr<CommonData> commonDataPtr;
206};
207
208template <int DIM> struct OpCalculateLogC_dC : public DomainEleOp {
209
211 boost::shared_ptr<CommonData> common_data)
213 commonDataPtr(common_data) {
214 std::fill(&doEntities[MBEDGE], &doEntities[MBMAXTYPE], false);
215 }
216
219
224
225 // const size_t nb_gauss_pts = matGradPtr->size2();
226 const size_t nb_gauss_pts = getGaussPts().size2();
227 constexpr auto size_symm = (DIM * (DIM + 1)) / 2;
228 commonDataPtr->matLogCdC.resize(size_symm * size_symm, nb_gauss_pts, false);
229 auto t_logC_dC = getFTensor4DdgFromMat<DIM, DIM>(commonDataPtr->matLogCdC);
230 auto t_eig_val = getFTensor1FromMat<DIM>(commonDataPtr->matEigVal);
231 auto t_eig_vec = getFTensor2FromMat<DIM, DIM>(commonDataPtr->matEigVec);
232
233 for (size_t gg = 0; gg != nb_gauss_pts; ++gg) {
234
237 eig(i) = t_eig_val(i);
238 eigen_vec(i, j) = t_eig_vec(i, j);
239
240 // rare case when two eigen values are equal
241 auto nb_uniq = get_uniq_nb<DIM>(&eig(0));
242 auto dlogC_dC = EigenMatrix::getDiffMat(eig, eigen_vec, f, d_f, nb_uniq);
243 dlogC_dC(i, j, k, l) *= 2;
244
245 t_logC_dC(i, j, k, l) = dlogC_dC(i, j, k, l);
246
247 ++t_logC_dC;
248 ++t_eig_val;
249 ++t_eig_vec;
250 }
251
253 }
254
255private:
256 boost::shared_ptr<CommonData> commonDataPtr;
257};
258
259template <int DIM> struct OpCalculateHenckyStress : public DomainEleOp {
260
262 boost::shared_ptr<CommonData> common_data)
264 commonDataPtr(common_data) {
265 std::fill(&doEntities[MBEDGE], &doEntities[MBMAXTYPE], false);
266 }
267
270
275
276 constexpr auto t_kd = FTensor::Kronecker_Delta<int>();
277
278 // const size_t nb_gauss_pts = matGradPtr->size2();
279 const size_t nb_gauss_pts = getGaussPts().size2();
280 auto t_D = getFTensor4DdgFromMat<DIM, DIM, 0>(*commonDataPtr->matDPtr);
281 auto t_logC = getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matLogC);
282 constexpr auto size_symm = (DIM * (DIM + 1)) / 2;
283 commonDataPtr->matHenckyStress.resize(size_symm, nb_gauss_pts, false);
284 auto t_T = getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matHenckyStress);
285
286 for (size_t gg = 0; gg != nb_gauss_pts; ++gg) {
287 t_T(i, j) = t_D(i, j, k, l) * t_logC(k, l);
288 ++t_logC;
289 ++t_T;
290 ++t_D;
291 }
292
294 }
295
296private:
297 boost::shared_ptr<CommonData> commonDataPtr;
298};
299
300template <int DIM> struct OpCalculateHenckyPlasticStress : public DomainEleOp {
301
303 boost::shared_ptr<CommonData> common_data,
304 boost::shared_ptr<MatrixDouble> mat_D_ptr,
305 const double scale = 1)
307 scaleStress(scale), matDPtr(mat_D_ptr) {
308 std::fill(&doEntities[MBEDGE], &doEntities[MBMAXTYPE], false);
309
310 matLogCPlastic = commonDataPtr->matLogCPlastic;
311 }
312
315
320
321 constexpr auto t_kd = FTensor::Kronecker_Delta<int>();
322
323 // const size_t nb_gauss_pts = matGradPtr->size2();
324 const size_t nb_gauss_pts = getGaussPts().size2();
325 auto t_D = getFTensor4DdgFromMat<DIM, DIM, 0>(*matDPtr);
326 auto t_logC = getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matLogC);
327 auto t_logCPlastic = getFTensor2SymmetricFromMat<DIM>(*matLogCPlastic);
328 constexpr auto size_symm = (DIM * (DIM + 1)) / 2;
329 commonDataPtr->matHenckyStress.resize(size_symm, nb_gauss_pts, false);
330 auto t_T = getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matHenckyStress);
331
332 for (size_t gg = 0; gg != nb_gauss_pts; ++gg) {
333 t_T(i, j) = t_D(i, j, k, l) * (t_logC(k, l) - t_logCPlastic(k, l));
334 t_T(i, j) /= scaleStress;
335 ++t_logC;
336 ++t_T;
337 ++t_D;
338 ++t_logCPlastic;
339 }
340
342 }
343
344private:
345 boost::shared_ptr<CommonData> commonDataPtr;
346 boost::shared_ptr<MatrixDouble> matDPtr;
347 boost::shared_ptr<MatrixDouble> matLogCPlastic;
348 const double scaleStress;
349};
350
351template <int DIM> struct OpCalculatePiolaStress : public DomainEleOp {
352
354 boost::shared_ptr<CommonData> common_data)
356 commonDataPtr(common_data) {
357 std::fill(&doEntities[MBEDGE], &doEntities[MBMAXTYPE], false);
358 }
359
362
367
368 constexpr auto t_kd = FTensor::Kronecker_Delta<int>();
369
370 // const size_t nb_gauss_pts = matGradPtr->size2();
371 const size_t nb_gauss_pts = getGaussPts().size2();
372 auto t_D = getFTensor4DdgFromMat<DIM, DIM, 0>(*commonDataPtr->matDPtr);
373 auto t_logC = getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matLogC);
374 auto t_logC_dC = getFTensor4DdgFromMat<DIM, DIM>(commonDataPtr->matLogCdC);
375 constexpr auto size_symm = (DIM * (DIM + 1)) / 2;
376 commonDataPtr->matFirstPiolaStress.resize(DIM * DIM, nb_gauss_pts, false);
377 commonDataPtr->matSecondPiolaStress.resize(size_symm, nb_gauss_pts, false);
378 auto t_P = getFTensor2FromMat<DIM, DIM>(commonDataPtr->matFirstPiolaStress);
379 auto t_T = getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matHenckyStress);
380 auto t_S =
381 getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matSecondPiolaStress);
382 auto t_grad = getFTensor2FromMat<DIM, DIM>(*(commonDataPtr->matGradPtr));
383
384 for (size_t gg = 0; gg != nb_gauss_pts; ++gg) {
386 t_F(i, j) = t_grad(i, j) + t_kd(i, j);
387 t_S(k, l) = t_T(i, j) * t_logC_dC(i, j, k, l);
388 t_P(i, l) = t_F(i, k) * t_S(k, l);
389
390 ++t_grad;
391 ++t_logC;
392 ++t_logC_dC;
393 ++t_P;
394 ++t_T;
395 ++t_S;
396 ++t_D;
397 }
398
400 }
401
402private:
403 boost::shared_ptr<CommonData> commonDataPtr;
404};
405
406template <int DIM> struct OpHenckyTangent : public DomainEleOp {
407 OpHenckyTangent(const std::string field_name,
408 boost::shared_ptr<CommonData> common_data,
409 boost::shared_ptr<MatrixDouble> mat_D_ptr = nullptr)
411 commonDataPtr(common_data) {
412 std::fill(&doEntities[MBEDGE], &doEntities[MBMAXTYPE], false);
413 if (mat_D_ptr)
414 matDPtr = mat_D_ptr;
415 else
416 matDPtr = commonDataPtr->matDPtr;
417 }
418
421
430
431 constexpr auto t_kd = FTensor::Kronecker_Delta<int>();
432 // const size_t nb_gauss_pts = matGradPtr->size2();
433 const size_t nb_gauss_pts = getGaussPts().size2();
434 constexpr auto size_symm = (DIM * (DIM + 1)) / 2;
435 commonDataPtr->matTangent.resize(DIM * DIM * DIM * DIM, nb_gauss_pts);
436 auto dP_dF =
437 getFTensor4FromMat<DIM, DIM, DIM, DIM, 1>(commonDataPtr->matTangent);
438
439 auto t_D = getFTensor4DdgFromMat<DIM, DIM, 0>(*matDPtr);
440 auto t_eig_val = getFTensor1FromMat<DIM>(commonDataPtr->matEigVal);
441 auto t_eig_vec = getFTensor2FromMat<DIM, DIM>(commonDataPtr->matEigVec);
442 auto t_T = getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matHenckyStress);
443 auto t_S =
444 getFTensor2SymmetricFromMat<DIM>(commonDataPtr->matSecondPiolaStress);
445 auto t_grad = getFTensor2FromMat<DIM, DIM>(*(commonDataPtr->matGradPtr));
446 auto t_logC_dC = getFTensor4DdgFromMat<DIM, DIM>(commonDataPtr->matLogCdC);
447
448 for (size_t gg = 0; gg != nb_gauss_pts; ++gg) {
449
451 t_F(i, j) = t_grad(i, j) + t_kd(i, j);
452
456 eig(i) = t_eig_val(i);
457 eigen_vec(i, j) = t_eig_vec(i, j);
458 T(i, j) = t_T(i, j);
459
460 // rare case when two eigen values are equal
461 auto nb_uniq = get_uniq_nb<DIM>(&eig(0));
462
464 dC_dF(i, j, k, l) = (t_kd(i, l) * t_F(k, j)) + (t_kd(j, l) * t_F(k, i));
465
466 auto TL =
467 EigenMatrix::getDiffDiffMat(eig, eigen_vec, f, d_f, dd_f, T, nb_uniq);
468
469 TL(i, j, k, l) *= 4;
470 FTensor::Ddg<double, DIM, DIM> P_D_P_plus_TL;
471 P_D_P_plus_TL(i, j, k, l) =
472 TL(i, j, k, l) +
473 (t_logC_dC(i, j, o, p) * t_D(o, p, m, n)) * t_logC_dC(m, n, k, l);
474 P_D_P_plus_TL(i, j, k, l) *= 0.5;
475 dP_dF(i, j, m, n) = t_kd(i, m) * (t_kd(k, n) * t_S(k, j));
476 dP_dF(i, j, m, n) +=
477 t_F(i, k) * (P_D_P_plus_TL(k, j, o, p) * dC_dF(o, p, m, n));
478
479 ++dP_dF;
480
481 ++t_grad;
482 ++t_eig_val;
483 ++t_eig_vec;
484 ++t_logC_dC;
485 ++t_S;
486 ++t_T;
487 ++t_D;
488 }
489
491 }
492
493private:
494 boost::shared_ptr<CommonData> commonDataPtr;
495 boost::shared_ptr<MatrixDouble> matDPtr;
496};
497
498template <int DIM> struct OpPostProcHencky : public DomainEleOp {
499 OpPostProcHencky(const std::string field_name,
500 moab::Interface &post_proc_mesh,
501 std::vector<EntityHandle> &map_gauss_pts,
502 boost::shared_ptr<CommonData> common_data_ptr);
503 MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
504
505private:
507 std::vector<EntityHandle> &mapGaussPts;
508 boost::shared_ptr<CommonData> commonDataPtr;
509};
510
511template <int DIM>
513 const std::string field_name, moab::Interface &post_proc_mesh,
514 std::vector<EntityHandle> &map_gauss_pts,
515 boost::shared_ptr<CommonData> common_data_ptr)
516 : DomainEleOp(field_name, DomainEleOp::OPROW), postProcMesh(post_proc_mesh),
517 mapGaussPts(map_gauss_pts), commonDataPtr(common_data_ptr) {
518 // Opetor is only executed for vertices
519 std::fill(&doEntities[MBEDGE], &doEntities[MBMAXTYPE], false);
520}
521
522//! [Postprocessing]
523template <int DIM>
525 EntData &data) {
527
531
532 std::array<double, 9> def;
533 std::fill(def.begin(), def.end(), 0);
534
535 auto get_tag = [&](const std::string name, size_t size) {
536 Tag th;
537 CHKERR postProcMesh.tag_get_handle(name.c_str(), size, MB_TYPE_DOUBLE, th,
538 MB_TAG_CREAT | MB_TAG_SPARSE,
539 def.data());
540 return th;
541 };
542
543 MatrixDouble3by3 mat(3, 3);
544
545 auto set_matrix = [&](auto &t) -> MatrixDouble3by3 & {
546 mat.clear();
547 for (size_t r = 0; r != SPACE_DIM; ++r)
548 for (size_t c = 0; c != SPACE_DIM; ++c)
549 mat(r, c) = t(r, c);
550 return mat;
551 };
552
553 auto set_tag = [&](auto th, auto gg, MatrixDouble3by3 &mat) {
554 return postProcMesh.tag_set_data(th, &mapGaussPts[gg], 1,
555 &*mat.data().begin());
556 };
557
558 auto th_grad = get_tag("GRAD", 9);
559 auto th_stress = get_tag("STRESS", 9);
560
561 auto t_piola =
562 getFTensor2FromMat<DIM, DIM>(commonDataPtr->matFirstPiolaStress);
563 auto t_grad = getFTensor2FromMat<DIM, DIM>(*commonDataPtr->matGradPtr);
564
567
568 for (int gg = 0; gg != commonDataPtr->matGradPtr->size2(); ++gg) {
569
570 F(i, j) = t_grad(i, j) + kronecker_delta(i, j);
571 double inv_t_detF;
572
573 if (DIM == 2)
574 determinantTensor2by2(F, inv_t_detF);
575 else
576 determinantTensor3by3(F, inv_t_detF);
577
578 inv_t_detF = 1. / inv_t_detF;
579
580 cauchy_stress(i, j) = inv_t_detF * (t_piola(i, k) ^ F(j, k));
581
582 CHKERR set_tag(th_grad, gg, set_matrix(t_grad));
583 CHKERR set_tag(th_stress, gg, set_matrix(cauchy_stress));
584
585 ++t_piola;
586 ++t_grad;
587 }
588
590}
591
592} // namespace HenckyOps
static Index< 'o', 3 > o
static Index< 'p', 3 > p
constexpr double a
constexpr int SPACE_DIM
Kronecker Delta class.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
FTensor::Index< 'n', SPACE_DIM > n
FTensor::Index< 'm', SPACE_DIM > m
constexpr auto t_kd
FTensor::Index< 'i', SPACE_DIM > i
const double v
phase velocity of light in medium (cm/ns)
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
const double T
FTensor::Ddg< double, 3, 3 > getDiffMat(Val< double, 3 > &t_val, Vec< double, 3 > &t_vec, Fun< double > f, Fun< double > d_f, const int nb)
Get the Diff Mat object.
FTensor::Tensor2_symmetric< double, 3 > getMat(Val< double, 3 > &t_val, Vec< double, 3 > &t_vec, Fun< double > f)
Get the Mat object.
FTensor::Ddg< double, 3, 3 > getDiffDiffMat(Val< double, 3 > &t_val, Vec< double, 3 > &t_vec, Fun< double > f, Fun< double > d_f, Fun< double > dd_f, FTensor::Tensor2< double, 3, 3 > &t_S, const int nb)
Tensor2_Expr< Kronecker_Delta< T >, T, Dim0, Dim1, i, j > kronecker_delta(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
Rank 2.
auto get_uniq_nb(double *ptr)
Definition: HenckyOps.hpp:15
auto sort_eigen_vals(FTensor::Tensor1< double, DIM > &eig, FTensor::Tensor2< double, DIM, DIM > &eigen_vec)
Definition: HenckyOps.hpp:23
auto f
Definition: HenckyOps.hpp:5
constexpr double eps
Definition: HenckyOps.hpp:3
auto dd_f
Definition: HenckyOps.hpp:7
auto d_f
Definition: HenckyOps.hpp:6
auto is_eq(const double &a, const double &b)
Definition: HenckyOps.hpp:11
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:77
MatrixBoundedArray< double, 9 > MatrixDouble3by3
Definition: Types.hpp:105
MoFEMErrorCode determinantTensor2by2(T1 &t, T2 &det)
Calculate determinant 2 by 2.
Definition: Templates.hpp:1240
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1955
static auto determinantTensor3by3(T &t)
Calculate the determinant of a 3x3 matrix or a tensor of rank 2.
Definition: Templates.hpp:1199
const double r
rate factor
double scale
Definition: plastic.cpp:94
constexpr double t
plate stiffness
Definition: plate.cpp:59
constexpr auto field_name
auto getMatFirstPiolaStress()
Definition: HenckyOps.hpp:76
MatrixDouble matEigVal
Definition: HenckyOps.hpp:67
MatrixDouble matLogCdC
Definition: HenckyOps.hpp:70
MatrixDouble matEigVec
Definition: HenckyOps.hpp:68
MatrixDouble matTangent
Definition: HenckyOps.hpp:74
boost::shared_ptr< MatrixDouble > matLogCPlastic
Definition: HenckyOps.hpp:65
MatrixDouble matFirstPiolaStress
Definition: HenckyOps.hpp:71
MatrixDouble matLogC
Definition: HenckyOps.hpp:69
boost::shared_ptr< MatrixDouble > matDPtr
Definition: HenckyOps.hpp:64
MatrixDouble matSecondPiolaStress
Definition: HenckyOps.hpp:72
MatrixDouble matHenckyStress
Definition: HenckyOps.hpp:73
boost::shared_ptr< MatrixDouble > matGradPtr
Definition: HenckyOps.hpp:63
boost::shared_ptr< CommonData > commonDataPtr
Definition: HenckyOps.hpp:155
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: HenckyOps.hpp:104
OpCalculateEigenVals(const std::string field_name, boost::shared_ptr< CommonData > common_data)
Definition: HenckyOps.hpp:97
OpCalculateHenckyPlasticStress(const std::string field_name, boost::shared_ptr< CommonData > common_data, boost::shared_ptr< MatrixDouble > mat_D_ptr, const double scale=1)
Definition: HenckyOps.hpp:302
boost::shared_ptr< MatrixDouble > matDPtr
Definition: HenckyOps.hpp:346
boost::shared_ptr< MatrixDouble > matLogCPlastic
Definition: HenckyOps.hpp:347
boost::shared_ptr< CommonData > commonDataPtr
Definition: HenckyOps.hpp:345
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: HenckyOps.hpp:313
OpCalculateHenckyStress(const std::string field_name, boost::shared_ptr< CommonData > common_data)
Definition: HenckyOps.hpp:261
boost::shared_ptr< CommonData > commonDataPtr
Definition: HenckyOps.hpp:297
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: HenckyOps.hpp:268
OpCalculateLogC_dC(const std::string field_name, boost::shared_ptr< CommonData > common_data)
Definition: HenckyOps.hpp:210
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: HenckyOps.hpp:217
boost::shared_ptr< CommonData > commonDataPtr
Definition: HenckyOps.hpp:256
boost::shared_ptr< CommonData > commonDataPtr
Definition: HenckyOps.hpp:205
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: HenckyOps.hpp:167
OpCalculateLogC(const std::string field_name, boost::shared_ptr< CommonData > common_data)
Definition: HenckyOps.hpp:160
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: HenckyOps.hpp:360
boost::shared_ptr< CommonData > commonDataPtr
Definition: HenckyOps.hpp:403
OpCalculatePiolaStress(const std::string field_name, boost::shared_ptr< CommonData > common_data)
Definition: HenckyOps.hpp:353
boost::shared_ptr< CommonData > commonDataPtr
Definition: HenckyOps.hpp:494
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: HenckyOps.hpp:419
boost::shared_ptr< MatrixDouble > matDPtr
Definition: HenckyOps.hpp:495
OpHenckyTangent(const std::string field_name, boost::shared_ptr< CommonData > common_data, boost::shared_ptr< MatrixDouble > mat_D_ptr=nullptr)
Definition: HenckyOps.hpp:407
Definition: HenckyOps.hpp:498
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
[Postprocessing]
Definition: HenckyOps.hpp:524
OpPostProcHencky(const std::string field_name, moab::Interface &post_proc_mesh, std::vector< EntityHandle > &map_gauss_pts, boost::shared_ptr< CommonData > common_data_ptr)
Definition: HenckyOps.hpp:512
std::vector< EntityHandle > & mapGaussPts
Definition: HenckyOps.hpp:507
moab::Interface & postProcMesh
Definition: HenckyOps.hpp:506
boost::shared_ptr< CommonData > commonDataPtr
Definition: HenckyOps.hpp:508
std::array< bool, MBMAXTYPE > doEntities
If true operator is executed for entity.
Data on single entity (This is passed as argument to DataOperator::doWork)
@ OPROW
operator doWork function is executed on FE rows
MatrixDouble & getGaussPts()
matrix of integration (Gauss) points for Volume Element