v0.14.0
Loading...
Searching...
No Matches
HookeInternalStressElement.hpp
Go to the documentation of this file.
1/**
2 * \file HookeInternalStressElement.hpp
3 * \example HookeInternalStressElement.hpp
4 *
5 * \brief Operators and data structures for calculating internal stress
6 *
7 */
8
9/*
10 * This file is part of MoFEM.
11 * MoFEM is free software: you can redistribute it and/or modify it under
12 * the terms of the GNU Lesser General Public License as published by the
13 * Free Software Foundation, either version 3 of the License, or (at your
14 * option) any later version.
15 *
16 * MoFEM is distributed in the hope that it will be useful, but WITHOUT
17 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
19 * License for more details.
20 *
21 * You should have received a copy of the GNU Lesser General Public
22 * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
23
24#ifndef __HOOKE_INTERNAL_STRESS_ELEMENT_HPP__
25#define __HOOKE_INTERNAL_STRESS_ELEMENT_HPP__
26
28
30
31 struct DataAtIntegrationPts : public HookeElement::DataAtIntegrationPts {
32
33 boost::shared_ptr<MatrixDouble> internalStressMat;
34 boost::shared_ptr<MatrixDouble> actualStressMat;
35 boost::shared_ptr<MatrixDouble> deviatoricStressMat;
36 boost::shared_ptr<MatrixDouble> hydrostaticStressMat;
37
38 boost::shared_ptr<MatrixDouble> spatPosMat;
39 boost::shared_ptr<MatrixDouble> meshNodePosMat;
40
42
43 internalStressMat = boost::shared_ptr<MatrixDouble>(new MatrixDouble());
44 actualStressMat = boost::shared_ptr<MatrixDouble>(new MatrixDouble());
45 deviatoricStressMat = boost::shared_ptr<MatrixDouble>(new MatrixDouble());
47 boost::shared_ptr<MatrixDouble>(new MatrixDouble());
48
49 spatPosMat = boost::shared_ptr<MatrixDouble>(new MatrixDouble());
50 meshNodePosMat = boost::shared_ptr<MatrixDouble>(new MatrixDouble());
51 }
52 };
53
55
56 OpGetInternalStress(const std::string row_field,
57 const std::string col_field,
58 boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
59 moab::Interface &input_mesh, char *stress_tag_name)
60 : VolUserDataOperator(row_field, col_field, OPROW, true),
61 dataAtPts(data_at_pts), inputMesh(input_mesh),
62 stressTagName(stress_tag_name) {
63 std::fill(&doEntities[MBEDGE], &doEntities[MBMAXTYPE], false);
64 };
65
66 MoFEMErrorCode doWork(int row_side, EntityType row_type, EntData &row_data);
67
68 protected:
69 boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
70 moab::Interface &inputMesh;
72 };
73
75
76 OpInternalStrain_dx(const std::string row_field,
77 boost::shared_ptr<DataAtIntegrationPts> data_at_pts)
78 : OpAssemble(row_field, row_field, data_at_pts, OPROW, true),
79 dataAtPts(data_at_pts){};
80
82
83 protected:
84 boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
85 };
86
87 template <int S = 0>
89
90 typedef boost::function<
91
93
95
96 )
97
98 >
100
102 const std::string row_field, const std::string col_field,
103 boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
104 StrainFunction strain_fun)
105 : VolUserDataOperator(row_field, col_field, OPROW, true),
106 dataAtPts(data_at_pts), strainFun(strain_fun) {
107 std::fill(&doEntities[MBEDGE], &doEntities[MBMAXTYPE], false);
108 }
109
110 MoFEMErrorCode doWork(int row_side, EntityType row_type, EntData &row_data);
111
112 protected:
114 boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
115 };
116
118 boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
119 map<int, BlockData>
120 &blockSetsPtr; // FIXME: (works only with the first block)
121 moab::Interface &outputMesh;
122 bool isALE;
126
127 OpSaveStress(const string row_field, const string col_field,
128 boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
129 map<int, BlockData> &block_sets_ptr,
130 moab::Interface &output_mesh, double scale_factor,
131 bool save_mean = false, bool is_ale = false,
132 bool is_field_disp = true)
133 : VolUserDataOperator(row_field, col_field, OPROW, true),
134 dataAtPts(data_at_pts), blockSetsPtr(block_sets_ptr),
135 outputMesh(output_mesh), scaleFactor(scale_factor), isALE(is_ale),
136 isFieldDisp(is_field_disp), saveMean(save_mean){};
137
138 MoFEMErrorCode doWork(int row_side, EntityType row_type, EntData &row_data);
139 };
140
141 template <class ELEMENT>
142 struct OpPostProcHookeElement : public ELEMENT::UserDataOperator {
143 boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
144 map<int, BlockData>
145 &blockSetsPtr; // FIXME: (works only with the first block)
146 moab::Interface &postProcMesh;
147 std::vector<EntityHandle> &mapGaussPts;
148 bool isALE;
150
151 OpPostProcHookeElement(const string row_field,
152 boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
153 map<int, BlockData> &block_sets_ptr,
154 moab::Interface &post_proc_mesh,
155 std::vector<EntityHandle> &map_gauss_pts,
156 bool is_ale = false, bool is_field_disp = true)
157 : ELEMENT::UserDataOperator(row_field, UserDataOperator::OPROW),
158 dataAtPts(data_at_pts), blockSetsPtr(block_sets_ptr),
159 postProcMesh(post_proc_mesh), mapGaussPts(map_gauss_pts),
160 isALE(is_ale), isFieldDisp(is_field_disp) {}
161
162 MoFEMErrorCode doWork(int side, EntityType type,
164 };
165};
166
168 int row_side, EntityType row_type, EntData &row_data) {
170
171 const int nb_integration_pts = row_data.getN().size1();
172
173 const EntityHandle ent = getFEEntityHandle();
174
175 const int val_num = 9 * nb_integration_pts;
176 std::vector<double> def_vals(val_num, 0.0);
177 Tag th_internal_stress;
178 CHKERR inputMesh.tag_get_handle(
179 stressTagName, val_num, MB_TYPE_DOUBLE, th_internal_stress,
180 MB_TAG_CREAT | MB_TAG_SPARSE, &*def_vals.begin());
181
182 dataAtPts->internalStressMat->resize(9, nb_integration_pts, false);
183
184 CHKERR inputMesh.tag_get_data(
185 th_internal_stress, &ent, 1,
186 &*(dataAtPts->internalStressMat->data().begin()));
187
189}
190
193 FTensor::Index<'i', 3> i;
194 FTensor::Index<'j', 3> j;
195 FTensor::Index<'k', 3> k;
196 FTensor::Index<'l', 3> l;
198
199 auto get_tensor1 = [](VectorDouble &v, const int r) {
201 &v(r + 0), &v(r + 1), &v(r + 2));
202 };
203
204 const int nb_integration_pts = getGaussPts().size2();
205
206 auto t_internal_stress =
207 getFTensor2FromMat<3, 3>(*dataAtPts->internalStressMat);
208
209 // get element volume
210 double vol = getVolume();
211 auto t_w = getFTensor0IntegrationWeight();
212
213 nF.resize(nbRows, false);
214 nF.clear();
215
216 // get derivatives of base functions on rows
217 auto t_row_diff_base = row_data.getFTensor1DiffN<3>();
218 const int row_nb_base_fun = row_data.getN().size2();
219
220 for (size_t gg = 0; gg != nb_integration_pts; ++gg) {
221
222 // calculate scalar weight times element volume
223 double a = t_w * vol;
224 auto t_nf = get_tensor1(nF, 0);
225
226 int rr = 0;
227 for (; rr != nbRows / 3; ++rr) {
228 t_nf(i) += a * t_row_diff_base(j) * t_internal_stress(i, j);
229 ++t_row_diff_base;
230 ++t_nf;
231 }
232
233 for (; rr != row_nb_base_fun; ++rr)
234 ++t_row_diff_base;
235
236 ++t_w;
237 ++t_internal_stress;
238 }
239
241}
242
243template <int S>
246 int row_side, EntityType row_type, EntData &row_data) {
247 FTensor::Index<'i', 3> i;
248 FTensor::Index<'j', 3> j;
249 FTensor::Index<'k', 3> k;
250 FTensor::Index<'l', 3> l;
252
253 auto tensor_to_tensor = [](const auto &t1, auto &t2) {
254 t2(0, 0) = t1(0, 0);
255 t2(1, 1) = t1(1, 1);
256 t2(2, 2) = t1(2, 2);
257 t2(0, 1) = t2(1, 0) = t1(1, 0);
258 t2(0, 2) = t2(2, 0) = t1(2, 0);
259 t2(1, 2) = t2(2, 1) = t1(2, 1);
260 };
261
262 const int nb_integration_pts = getGaussPts().size2();
263 auto t_coords = getFTensor1CoordsAtGaussPts();
264
265 // elastic stiffness tensor (4th rank tensor with minor and major
266 // symmetry)
268 MAT_TO_DDG(dataAtPts->stiffnessMat));
269
270 dataAtPts->internalStressMat->resize(9, nb_integration_pts, false);
271 auto t_internal_stress =
272 getFTensor2FromMat<3, 3>(*dataAtPts->internalStressMat);
273
274 FTensor::Tensor2_symmetric<double, 3> t_internal_stress_symm;
275
276 for (size_t gg = 0; gg != nb_integration_pts; ++gg) {
277
278 auto t_fun_strain = strainFun(t_coords);
279 t_internal_stress_symm(i, j) = -t_D(i, j, k, l) * t_fun_strain(k, l);
280
281 tensor_to_tensor(t_internal_stress_symm, t_internal_stress);
282
283 ++t_coords;
284 ++t_D;
285 ++t_internal_stress;
286 }
287
289}
290
292 int row_side, EntityType row_type, EntData &row_data) {
294
295 if (row_type != MBVERTEX) {
297 }
298
299 const EntityHandle ent = getFEEntityHandle();
300
301 auto tensor_to_tensor = [](const auto &t1, auto &t2) {
302 t2(0, 0) = t1(0, 0);
303 t2(1, 1) = t1(1, 1);
304 t2(2, 2) = t1(2, 2);
305 t2(0, 1) = t2(1, 0) = t1(1, 0);
306 t2(0, 2) = t2(2, 0) = t1(2, 0);
307 t2(1, 2) = t2(2, 1) = t1(2, 1);
308 };
309
310 auto tensor_to_vector = [](const auto &t, auto &v) {
311 v(0) = t(0, 0);
312 v(1) = t(1, 1);
313 v(2) = t(2, 2);
314 v(3) = t(0, 1);
315 v(4) = t(0, 2);
316 v(5) = t(1, 2);
317 };
318
319 auto get_tag_handle = [&](auto name, auto size) {
320 Tag th;
321 std::vector<double> def_vals(size, 0.0);
322 CHKERR outputMesh.tag_get_handle(name, size, MB_TYPE_DOUBLE, th,
323 MB_TAG_CREAT | MB_TAG_SPARSE,
324 def_vals.data());
325 return th;
326 };
327
328 const int nb_integration_pts = row_data.getN().size1();
329
330 auto th_internal_stress =
331 get_tag_handle("INTERNAL_STRESS", 9 * nb_integration_pts);
332 auto th_actual_stress =
333 get_tag_handle("ACTUAL_STRESS", 9 * nb_integration_pts);
334 auto th_deviatoric_stress =
335 get_tag_handle("DEVIATORIC_STRESS", 9 * nb_integration_pts);
336 auto th_hydrostatic_stress =
337 get_tag_handle("HYDROSTATIC_STRESS", 9 * nb_integration_pts);
338
339 auto th_internal_stress_mean = get_tag_handle("MED_INTERNAL_STRESS", 9);
340 auto th_actual_stress_mean = get_tag_handle("MED_ACTUAL_STRESS", 9);
341
342 FTensor::Index<'i', 3> i;
343 FTensor::Index<'j', 3> j;
344 FTensor::Index<'k', 3> k;
345 FTensor::Index<'l', 3> l;
346
347 auto t_h = getFTensor2FromMat<3, 3>(*dataAtPts->hMat);
348 auto t_H = getFTensor2FromMat<3, 3>(*dataAtPts->HMat);
349
350 dataAtPts->stiffnessMat->resize(36, 1, false);
352 MAT_TO_DDG(dataAtPts->stiffnessMat));
353 for (auto &m : (blockSetsPtr)) {
354 const double young = m.second.E;
355 const double poisson = m.second.PoissonRatio;
356
357 const double coefficient = young / ((1 + poisson) * (1 - 2 * poisson));
358
359 t_D(i, j, k, l) = 0.;
360 t_D(0, 0, 0, 0) = t_D(1, 1, 1, 1) = t_D(2, 2, 2, 2) = 1 - poisson;
361 t_D(0, 1, 0, 1) = t_D(0, 2, 0, 2) = t_D(1, 2, 1, 2) =
362 0.5 * (1 - 2 * poisson);
363 t_D(0, 0, 1, 1) = t_D(1, 1, 0, 0) = t_D(0, 0, 2, 2) = t_D(2, 2, 0, 0) =
364 t_D(1, 1, 2, 2) = t_D(2, 2, 1, 1) = poisson;
365 t_D(i, j, k, l) *= coefficient;
366
367 break; // FIXME: calculates only first block
368 }
369
370 double detH = 0.;
373
375 FTensor::Tensor2_symmetric<double, 3> t_small_strain_symm;
376
378
379 auto t_internal_stress =
380 getFTensor2FromMat<3, 3>(*dataAtPts->internalStressMat);
381
382 dataAtPts->actualStressMat->resize(9, nb_integration_pts, false);
383 auto t_actual_stress = getFTensor2FromMat<3, 3>(*dataAtPts->actualStressMat);
384
385 dataAtPts->deviatoricStressMat->resize(9, nb_integration_pts, false);
386 auto t_deviatoric_stress =
387 getFTensor2FromMat<3, 3>(*dataAtPts->deviatoricStressMat);
388 dataAtPts->hydrostaticStressMat->resize(9, nb_integration_pts, false);
389 auto t_hydrostatic_stress =
390 getFTensor2FromMat<3, 3>(*dataAtPts->hydrostaticStressMat);
391
392 FTensor::Tensor2<double, 3, 3> t_internal_stress_mean;
393 FTensor::Tensor2<double, 3, 3> t_actual_stress_mean;
394
395 t_internal_stress_mean(i, j) = 0.;
396 t_actual_stress_mean(i, j) = 0.;
397
398 auto t_w = getFTensor0IntegrationWeight();
399
400 for (int gg = 0; gg != nb_integration_pts; ++gg) {
401
402 if (!isALE) {
403 t_small_strain_symm(i, j) = (t_h(i, j) || t_h(j, i)) / 2.;
404 } else {
405 CHKERR determinantTensor3by3(t_H, detH);
406 CHKERR invertTensor3by3(t_H, detH, t_invH);
407 t_F(i, j) = t_h(i, k) * t_invH(k, j);
408 t_small_strain_symm(i, j) = (t_F(i, j) || t_F(j, i)) / 2.;
409 ++t_H;
410 }
411
412 if (isALE || !isFieldDisp) {
413 for (auto ii : {0, 1, 2}) {
414 t_small_strain_symm(ii, ii) -= 1.;
415 }
416 }
417
418 // symmetric tensors need improvement
419 t_stress_symm(i, j) = t_D(i, j, k, l) * t_small_strain_symm(k, l);
420 tensor_to_tensor(t_stress_symm, t_stress);
421 t_actual_stress(i, j) = t_stress(i, j) + t_internal_stress(i, j);
422
423 t_actual_stress(i, j) *= scaleFactor;
424 t_internal_stress(i, j) *= scaleFactor;
425
426 double hydrostatic_pressure =
427 (t_actual_stress(0, 0) + t_actual_stress(1, 1) +
428 t_actual_stress(2, 2)) /
429 3.;
430
431 t_hydrostatic_stress(i, j) = 0.;
432 for (auto ii : {0, 1, 2}) {
433 t_hydrostatic_stress(ii, ii) += hydrostatic_pressure;
434 }
435
436 t_deviatoric_stress(i, j) = t_actual_stress(i, j);
437 for (auto ii : {0, 1, 2}) {
438 t_deviatoric_stress(ii, ii) -= hydrostatic_pressure;
439 }
440
441 t_actual_stress_mean(i, j) += t_w * t_actual_stress(i, j);
442 t_internal_stress_mean(i, j) += t_w * t_internal_stress(i, j);
443
444 ++t_w;
445 ++t_h;
446 ++t_actual_stress;
447 ++t_internal_stress;
448 }
449
450 if (saveMean) {
451 VectorDouble vec_actual_stress_mean;
452 vec_actual_stress_mean.resize(9, false);
453 vec_actual_stress_mean.clear();
454
455 VectorDouble vec_internal_stress_mean;
456 vec_internal_stress_mean.resize(9, false);
457 vec_internal_stress_mean.clear();
458
459 tensor_to_vector(t_actual_stress_mean, vec_actual_stress_mean);
460 tensor_to_vector(t_internal_stress_mean, vec_internal_stress_mean);
461
462 CHKERR outputMesh.tag_set_data(th_actual_stress_mean, &ent, 1,
463 &*(vec_actual_stress_mean.data().begin()));
464 CHKERR outputMesh.tag_set_data(th_internal_stress_mean, &ent, 1,
465 &*(vec_internal_stress_mean.data().begin()));
466 } else {
467 CHKERR outputMesh.tag_set_data(
468 th_internal_stress, &ent, 1,
469 &*(dataAtPts->internalStressMat->data().begin()));
470 CHKERR outputMesh.tag_set_data(
471 th_actual_stress, &ent, 1,
472 &*(dataAtPts->actualStressMat->data().begin()));
473
474 CHKERR outputMesh.tag_set_data(
475 th_hydrostatic_stress, &ent, 1,
476 &*(dataAtPts->hydrostaticStressMat->data().begin()));
477 CHKERR outputMesh.tag_set_data(
478 th_deviatoric_stress, &ent, 1,
479 &*(dataAtPts->deviatoricStressMat->data().begin()));
480 }
481
483}
484
485template <class ELEMENT>
488 int side, EntityType type, EntitiesFieldData::EntData &data) {
490
491 if (type != MBVERTEX) {
493 }
494
495 auto tensor_to_tensor = [](const auto &t1, auto &t2) {
496 t2(0, 0) = t1(0, 0);
497 t2(1, 1) = t1(1, 1);
498 t2(2, 2) = t1(2, 2);
499 t2(0, 1) = t2(1, 0) = t1(1, 0);
500 t2(0, 2) = t2(2, 0) = t1(2, 0);
501 t2(1, 2) = t2(2, 1) = t1(2, 1);
502 };
503
504 auto get_tag_handle = [&](auto name, auto size) {
505 Tag th;
506 std::vector<double> def_vals(size, 0.0);
507 CHKERR postProcMesh.tag_get_handle(name, size, MB_TYPE_DOUBLE, th,
508 MB_TAG_CREAT | MB_TAG_SPARSE,
509 def_vals.data());
510 return th;
511 };
512
513 auto th_stress = get_tag_handle("STRESS", 9);
514 auto th_strain = get_tag_handle("STRAIN", 9);
515 auto th_psi = get_tag_handle("ENERGY", 1);
516 auto th_disp = get_tag_handle("DISPLACEMENT", 3);
517
518 const int nb_integration_pts = mapGaussPts.size();
519
520 FTensor::Index<'i', 3> i;
521 FTensor::Index<'j', 3> j;
522 FTensor::Index<'k', 3> k;
523 FTensor::Index<'l', 3> l;
524
525 auto t_h = getFTensor2FromMat<3, 3>(*dataAtPts->hMat);
526 auto t_H = getFTensor2FromMat<3, 3>(*dataAtPts->HMat);
527
528 dataAtPts->stiffnessMat->resize(36, 1, false);
530 MAT_TO_DDG(dataAtPts->stiffnessMat));
531 for (auto &m : (blockSetsPtr)) {
532 const double young = m.second.E;
533 const double poisson = m.second.PoissonRatio;
534
535 const double coefficient = young / ((1 + poisson) * (1 - 2 * poisson));
536
537 t_D(i, j, k, l) = 0.;
538 t_D(0, 0, 0, 0) = t_D(1, 1, 1, 1) = t_D(2, 2, 2, 2) = 1 - poisson;
539 t_D(0, 1, 0, 1) = t_D(0, 2, 0, 2) = t_D(1, 2, 1, 2) =
540 0.5 * (1 - 2 * poisson);
541 t_D(0, 0, 1, 1) = t_D(1, 1, 0, 0) = t_D(0, 0, 2, 2) = t_D(2, 2, 0, 0) =
542 t_D(1, 1, 2, 2) = t_D(2, 2, 1, 1) = poisson;
543 t_D(i, j, k, l) *= coefficient;
544
545 break; // FIXME: calculates only first block
546 }
547
548 double detH = 0.;
552 FTensor::Tensor2<double, 3, 3> t_small_strain;
553
555 FTensor::Tensor2_symmetric<double, 3> t_small_strain_symm;
556
557 auto t_spat_pos = getFTensor1FromMat<3>(*dataAtPts->spatPosMat);
558 auto t_mesh_node_pos = getFTensor1FromMat<3>(*dataAtPts->meshNodePosMat);
560
561 for (int gg = 0; gg != nb_integration_pts; ++gg) {
562
563 if (isFieldDisp) {
564 t_h(0, 0) += 1;
565 t_h(1, 1) += 1;
566 t_h(2, 2) += 1;
567 }
568
569 if (!isALE) {
570 t_small_strain_symm(i, j) = (t_h(i, j) || t_h(j, i)) / 2.;
571 } else {
572 CHKERR determinantTensor3by3(t_H, detH);
573 CHKERR invertTensor3by3(t_H, detH, t_invH);
574 t_F(i, j) = t_h(i, k) * t_invH(k, j);
575 t_small_strain_symm(i, j) = (t_F(i, j) || t_F(j, i)) / 2.;
576 ++t_H;
577 }
578
579 t_small_strain_symm(0, 0) -= 1;
580 t_small_strain_symm(1, 1) -= 1;
581 t_small_strain_symm(2, 2) -= 1;
582
583 // symmetric tensors need improvement
584 t_stress_symm(i, j) = t_D(i, j, k, l) * t_small_strain_symm(k, l);
585 tensor_to_tensor(t_stress_symm, t_stress);
586 tensor_to_tensor(t_small_strain_symm, t_small_strain);
587
588 t_disp(i) = t_spat_pos(i) - t_mesh_node_pos(i);
589
590 const double psi = 0.5 * t_stress_symm(i, j) * t_small_strain_symm(i, j);
591
592 CHKERR postProcMesh.tag_set_data(th_psi, &mapGaussPts[gg], 1, &psi);
593 CHKERR postProcMesh.tag_set_data(th_stress, &mapGaussPts[gg], 1,
594 &t_stress(0, 0));
595 CHKERR postProcMesh.tag_set_data(th_strain, &mapGaussPts[gg], 1,
596 &t_small_strain(0, 0));
597 CHKERR postProcMesh.tag_set_data(th_disp, &mapGaussPts[gg], 1, &t_disp(0));
598
599 ++t_h;
600 ++t_spat_pos;
601 ++t_mesh_node_pos;
602 }
603
605}
606
607#endif
#define MAT_TO_DDG(SM)
ForcesAndSourcesCore::UserDataOperator UserDataOperator
constexpr double a
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
FTensor::Index< 'm', SPACE_DIM > m
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
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
UBlasMatrix< double > MatrixDouble
Definition Types.hpp:77
MoFEMErrorCode invertTensor3by3(ublas::matrix< T, L, A > &jac_data, ublas::vector< T, A > &det_data, ublas::matrix< T, L, A > &inv_jac_data)
Calculate inverse of tensor rank 2 at integration points.
static auto determinantTensor3by3(T &t)
Calculate the determinant of a 3x3 matrix or a tensor of rank 2.
constexpr double t
plate stiffness
Definition plate.cpp:59
OpGetAnalyticalInternalStress(const std::string row_field, const std::string col_field, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, StrainFunction strain_fun)
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntData &row_data)
Operator for linear form, usually to calculate values on right hand side.
boost::function< FTensor::Tensor2_symmetric< double, 3 >(FTensor::Tensor1< FTensor::PackPtr< double *, 3 >, 3 > &t_coords) StrainFunction)
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
OpGetInternalStress(const std::string row_field, const std::string col_field, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, moab::Interface &input_mesh, char *stress_tag_name)
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntData &row_data)
Operator for linear form, usually to calculate values on right hand side.
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
OpInternalStrain_dx(const std::string row_field, boost::shared_ptr< DataAtIntegrationPts > data_at_pts)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
bool isALE
bool isFieldDisp
map< int, BlockData > & blockSetsPtr
moab::Interface & postProcMesh
OpPostProcHookeElement(const string row_field, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, map< int, BlockData > &block_sets_ptr, moab::Interface &post_proc_mesh, std::vector< EntityHandle > &map_gauss_pts, bool is_ale=false, bool is_field_disp=true)
std::vector< EntityHandle > & mapGaussPts
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
OpSaveStress(const string row_field, const string col_field, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, map< int, BlockData > &block_sets_ptr, moab::Interface &output_mesh, double scale_factor, bool save_mean=false, bool is_ale=false, bool is_field_disp=true)
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntData &row_data)
Operator for linear form, usually to calculate values on right hand side.
std::array< bool, MBMAXTYPE > doEntities
If true operator is executed for entity.
Data on single entity (This is passed as argument to DataOperator::doWork)
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1DiffN(const FieldApproximationBase base)
Get derivatives of base functions.
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
EntityHandle getFEEntityHandle() const
Return finite element entity handle.
@ OPROW
operator doWork function is executed on FE rows