v0.14.0
Loading...
Searching...
No Matches
RotatingFrameOperators.cpp
Go to the documentation of this file.
1/** \file RotatingFrameOperators.cpp
2 */
3
4/* This file is part of MoFEM.
5 * MoFEM is free software: you can redistribute it and/or modify it under
6 * the terms of the GNU Lesser General Public License as published by the
7 * Free Software Foundation, either version 3 of the License, or (at your
8 * option) any later version.
9 *
10 * MoFEM is distributed in the hope that it will be useful, but WITHOUT
11 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
12 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
13 * License for more details.
14 *
15 * You should have received a copy of the GNU Lesser General Public
16 * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
17
18#include <MoFEM.hpp>
19using namespace MoFEM;
20using namespace FTensor;
21
24#include <RigidBodies.hpp>
25#include <BlockMatData.hpp>
27#include <BasicFeTools.hpp>
28#include <MatrixFunction.hpp>
29
30#include <ElasticOperators.hpp>
31#include <PlasticOperators.hpp>
33
34#include <chrono> //for debugging
35
36struct MeasureTime {
37 int time;
38 chrono::high_resolution_clock::time_point start;
39 chrono::high_resolution_clock::time_point stop;
40 MeasureTime() { start = chrono::high_resolution_clock::now(); }
42 stop = chrono::high_resolution_clock::now();
43 auto duration = chrono::duration_cast<chrono::microseconds>(stop - start);
44 cout << "Time taken by function: " << duration.count() << " ms." << endl;
45 }
46};
47
49
51 const std::string field_name, boost::shared_ptr<CommonData> common_data_ptr)
53 commonDataPtr(common_data_ptr) {}
54
57
58 const size_t nb_dofs = data.getIndices().size();
59
60 if (nb_dofs) {
61
62 const size_t nb_base_functions = data.getN().size2();
63 if (3 * nb_base_functions < nb_dofs)
64 SETERRQ(
65 PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
66 "Number of DOFs is larger than number of base functions on entity");
67
68 const size_t nb_gauss_pts = data.getN().size1();
69
70 auto t_w = getFTensor0IntegrationWeight();
71 auto t_diff_base = data.getFTensor1DiffN<3>();
72 auto t_grad = getFTensor2FromMat<3, 3>(*(commonDataPtr->mGradPtr));
73 auto t_coords = getFTensor1CoordsAtGaussPts();
74
75 for (size_t gg = 0; gg != nb_gauss_pts; ++gg) {
76
78 omg(i) = (*cache).Omega(i, j) * t_coords(j);
79
80 double alpha = getMeasure() * t_w * (*cache).density;
81 auto t_nf = getFTensor1FromArray<3, 3>(locF);
82
83 size_t bb = 0;
84 for (; bb != nb_dofs / 3; ++bb) {
85
86 t_nf(i) -= alpha * (t_grad(i, j) + kronecker_delta(i, j)) * omg(j) *
87 t_diff_base(k) * omg(k);
88
89 ++t_diff_base;
90 ++t_nf;
91 }
92 for (; bb < nb_base_functions; ++bb)
93 ++t_diff_base;
94
95 ++t_coords;
96 ++t_grad;
97 ++t_w;
98 }
99 }
100
102}
103
105 const std::string row_field_name, const std::string col_field_name,
106 boost::shared_ptr<CommonData> common_data_ptr)
107 : DomainEleOpAssembly(row_field_name, col_field_name,
108 DomainEleOp::OPROWCOL),
109 commonDataPtr(common_data_ptr) {
110 sYmm = false;
111}
112
113//! [Stiffness]
115 EntData &col_data) {
117
118 const size_t nb_row_dofs = row_data.getIndices().size();
119 const size_t nb_col_dofs = col_data.getIndices().size();
120
121 if (nb_row_dofs && nb_col_dofs) {
122
123 const size_t nb_integration_pts = row_data.getN().size1();
124 const size_t nb_row_base_funcs = row_data.getN().size2();
125 auto t_row_diff_base = row_data.getFTensor1DiffN<3>();
126 auto t_w = getFTensor0IntegrationWeight();
127 auto t_coords = getFTensor1CoordsAtGaussPts();
128 for (size_t gg = 0; gg != nb_integration_pts; ++gg) {
129 double alpha = getMeasure() * t_w * (*cache).density;
130
132 omg(i) = (*cache).Omega(i, j) * t_coords(j);
133
134 size_t rr = 0;
135 for (; rr != nb_row_dofs / 3; ++rr) {
136
137 auto t_mat = getFTensor2FromArray<3, 3, 3>(locMat, 3 * rr);
138
139 auto t_col_diff_base = col_data.getFTensor1DiffN<3>(gg, 0);
140
141 for (size_t cc = 0; cc != nb_col_dofs / 3; ++cc) {
142 t_mat(i, k) -= alpha * t_row_diff_base(j) * omg(j) *
143 t_col_diff_base(l) * omg(l) * kronecker_delta(i, k);
144
145 ++t_col_diff_base;
146 ++t_mat;
147 }
148
149 ++t_row_diff_base;
150 }
151 for (; rr != nb_row_base_funcs; ++rr)
152 ++t_row_diff_base;
153
154 ++t_coords;
155 ++t_w;
156 }
157 }
158
160}
161
163 const std::string field_name, boost::shared_ptr<CommonData> common_data_ptr,
164 boost::shared_ptr<DomainSideEle> side_op_fe)
166 commonDataPtr(common_data_ptr), sideOpFe(side_op_fe) {}
167
170
171 const size_t nb_gauss_pts = getGaussPts().size2();
172 const size_t nb_dofs = data.getIndices().size();
173
174 if (nb_dofs) {
175
176 auto t_w = getFTensor0IntegrationWeight();
177 auto t_coords = getFTensor1CoordsAtGaussPts();
178 CHKERR loopSideVolumes("dFE", *sideOpFe);
179 auto t_grad = getFTensor2FromMat<3, 3>(*(commonDataPtr->mGradPtr));
180 size_t nb_base_functions = data.getN().size2();
181 auto t_base = data.getFTensor0N();
182 auto t_normal = getFTensor1Normal();
183 const double ls = sqrt(t_normal(i) * t_normal(i));
184 t_normal(i) /= ls;
185
186 for (size_t gg = 0; gg != nb_gauss_pts; ++gg) {
187
188 auto t_nf = getFTensor1FromArray<3, 3>(locF);
189
190 if (getNormalsAtGaussPts().size1() == nb_gauss_pts) {
191 VectorDouble n = getNormalsAtGaussPts(gg);
192 auto t_n = getFTensor1FromPtr<3>(&*n.data().begin());
193 t_normal(i) = t_n(i) / sqrt(t_n(j) * t_n(j));
194 }
195
196 const double alpha = t_w * getMeasure() * (*cache).density;
198 omg(i) = (*cache).Omega(i, j) * t_coords(j);
199
200 size_t bb = 0;
201 for (; bb != nb_dofs / 3; ++bb) {
202 // cerr << test << endl;
203 t_nf(i) += alpha * t_base *
204 ((t_grad(i, j) + kronecker_delta(i, j)) * omg(j)) *
205 (omg(k) * t_normal(k));
206 // alpha * t_base * (t_grad(i, j) * omg(j));
207
208 ++t_nf;
209 ++t_base;
210 }
211 for (; bb < nb_base_functions; ++bb)
212 ++t_base;
213
214 ++t_grad;
215 ++t_coords;
216 ++t_w;
217 }
218 }
219
221}
222
224 const std::string row_field_name, const std::string col_field_name,
225 boost::shared_ptr<CommonData> common_data_ptr,
226 boost::shared_ptr<DomainSideEle> side_op_fe)
227 : BoundaryEleOp(row_field_name, col_field_name, OPROWCOL),
228 commonDataPtr(common_data_ptr), sideOpFe(side_op_fe) {
229 sYmm = false;
230}
231
233 EntityType row_type,
234 EntityType col_type,
235 EntData &row_data,
236 EntData &col_data) {
238
239 const size_t nb_gauss_pts = getGaussPts().size2();
240 const size_t row_nb_dofs = row_data.getIndices().size();
241 size_t col_nb_dofs = col_data.getIndices().size();
242
243 if (row_nb_dofs && col_nb_dofs) {
244
245 if (col_type == MBVERTEX) {
246 commonDataPtr->faceRowData = &row_data;
247 commonDataPtr->mMeasure = getMeasure();
249 }
250 }
251
253}
254
256 const std::string field_name, boost::shared_ptr<CommonData> common_data_ptr)
258 commonDataPtr(common_data_ptr) {
259 sYmm = false;
260}
261
263 EntData &col_data) {
264 return 0;
265}
266
268 int row_side, int col_side, EntityType row_type, EntityType col_type,
270 EntitiesFieldData::EntData &col_data) {
272 if (row_type != MBVERTEX)
274
275 const size_t nb_gauss_pts = getGaussPts().size2();
276 const size_t row_nb_dofs = commonDataPtr->faceRowData->getIndices().size();
277 const size_t col_nb_dofs = col_data.getIndices().size();
278
279 if (row_nb_dofs && col_nb_dofs) {
280
281 auto t_coords = getFTensor1CoordsAtGaussPts();
282 auto t_w = getFTensor0IntegrationWeight();
283 auto t_row_base = commonDataPtr->faceRowData->getFTensor0N();
284
285 auto t_row_diff_base = commonDataPtr->faceRowData->getFTensor1DiffN<3>();
286 size_t nb_face_functions = commonDataPtr->faceRowData->getN().size2();
287 auto t_normal = getFTensor1Normal();
288 const double ls = sqrt(t_normal(i) * t_normal(i));
289 t_normal(i) /= ls;
290
291 locMat.resize(row_nb_dofs, col_nb_dofs, false);
292 locMat.clear();
293
294 for (size_t gg = 0; gg != nb_gauss_pts; ++gg) {
295 // we are actually integrating on a triangle here, therefore 0.5
296 // const double alpha_test = t_w * getMeasure() * (*cache).density;
297 const double alpha = t_w * commonDataPtr->mMeasure * (*cache).density;
298 // const double alpha = t_w * 0.5 * (*cache).density;
300 omg(i) = (*cache).Omega(i, j) * t_coords(j);
301
302 if (getNormalsAtGaussPts().size1() == nb_gauss_pts) {
303 VectorDouble n = getNormalsAtGaussPts(gg);
304 auto t_n = getFTensor1FromPtr<3>(&*n.data().begin());
305 t_normal(i) = t_n(i) / sqrt(t_n(j) * t_n(j));
306 }
307
308 size_t rr = 0;
309 for (; rr != row_nb_dofs / 3; ++rr) {
310 auto t_mat = getFTensor2FromArray<3, 3, 3>(locMat, 3 * rr);
311
312 const double beta = alpha * t_row_base;
313 auto t_col_diff_base = col_data.getFTensor1DiffN<3>(gg, 0);
314
315 for (size_t cc = 0; cc != col_nb_dofs / 3; ++cc) {
316
317 t_mat(i, k) += beta * kronecker_delta(i, k) *
318 (omg(j) * t_col_diff_base(j)) * (omg(l) * t_normal(l));
319
320 ++t_col_diff_base;
321 ++t_mat;
322 }
323
324 ++t_row_base;
325 }
326 for (; rr < nb_face_functions; ++rr)
327 ++t_row_base;
328
329 ++t_coords;
330 ++t_w;
331 }
332 // CHKERR iNtegrate(row_data, col_data);
333 CHKERR aSsemble(*commonDataPtr->faceRowData, col_data, false);
334 }
336}
337
339 const std::string field_name, boost::shared_ptr<CommonData> common_data_ptr)
341 commonDataPtr(common_data_ptr) {}
342
344 EntityType type,
345 EntData &data) {
347 const size_t nb_dofs = data.getIndices().size();
348 if (nb_dofs) {
349 const size_t nb_integration_pts = data.getN().size1();
350
351 auto t_plastic_strain_dot =
352 getFTensor2SymmetricFromMat<3>(*(commonDataPtr->plasticStrainDotPtr));
353 auto t_tau_dot = getFTensor0FromVec(*(commonDataPtr->plasticTauDotPtr));
354 auto t_coords = getFTensor1CoordsAtGaussPts();
355 auto t_grad_tau =
356 getFTensor1FromMat<3>(*(commonDataPtr->plasticGradTauPtr));
357
358 // OpCalculateTensor2SymmetricFieldGradient
359 auto t_grad_plastic_strain =
360 getFTensor3DgFromMat<3, 3>(*(commonDataPtr->plasticGradStrainPtr));
361 commonDataPtr->guidingVelocityPtr->resize(3, nb_integration_pts, false);
362 commonDataPtr->guidingVelocityPtr->clear();
363 auto t_omega = getFTensor1FromMat<3>(*commonDataPtr->guidingVelocityPtr);
364
365 for (size_t gg = 0; gg != nb_integration_pts; ++gg) {
366
367 t_omega(i) = (*cache).Omega(i, j) * t_coords(j);
368 t_tau_dot += t_grad_tau(i) * t_omega(i);
369 t_plastic_strain_dot(i, j) += t_grad_plastic_strain(i, j, k) * t_omega(k);
370
371 ++t_coords;
372 ++t_omega;
373 ++t_grad_plastic_strain;
374 ++t_plastic_strain_dot;
375 ++t_tau_dot;
376 ++t_grad_tau;
377 }
378 }
379
381}
382
384 const std::string field_name, boost::shared_ptr<CommonData> common_data_ptr,
385 boost::shared_ptr<DomainSideEle> side_op_fe)
387 commonDataPtr(common_data_ptr), sideOpFe(side_op_fe) {
388 std::fill(&doEntities[MBVERTEX], &doEntities[MBMAXTYPE], false);
389 doEntities[MBTRI] = doEntities[MBQUAD] = true;
390}
391
393 EntData &data) {
395 // if (skipEnts.find(getFEEntityHandle()) != skipEnts.end())
396 // MoFEMFunctionReturnHot(0);
397 if ((type == MBTRI || type == MBQUAD) && side == 0) {
398 const size_t nb_dofs = data.getIndices().size();
399 // if (nb_dofs == 0)
400 // MoFEMFunctionReturnHot(0);
401
402 const size_t nb_integration_pts = data.getN().size1();
403 const size_t nb_base_functions = data.getN().size2();
404 commonDataPtr->plasticStrainSideMap.clear();
405 CHKERR loopSideVolumes("dFE", *sideOpFe);
406 const int check_size = commonDataPtr->plasticStrainSideMap.size();
407
408 if (commonDataPtr->plasticStrainSideMap.size() == 1)
409 if (static_cast<int>(
410 commonDataPtr->plasticStrainSideMap.begin()->second.size2()) !=
411 nb_integration_pts)
412 SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
413 "wrong number of integration points");
414
415 auto &plastic_mat_l = commonDataPtr->plasticStrainSideMap.at(-1);
416 auto &plastic_mat_r = commonDataPtr->plasticStrainSideMap.at(1);
417 auto t_plastic_strain_l = getFTensor2SymmetricFromMat<3>(plastic_mat_l);
418 auto t_plastic_strain_r = getFTensor2SymmetricFromMat<3>(plastic_mat_r);
419 auto &tau_vec_l = commonDataPtr->plasticTauSideMap.at(-1);
420 auto &tau_vec_r = commonDataPtr->plasticTauSideMap.at(1);
421 auto t_tau_l = getFTensor0FromVec(tau_vec_l);
422 auto t_tau_r = getFTensor0FromVec(tau_vec_r);
423
424 auto t_w = getFTensor0IntegrationWeight();
425 auto t_base = data.getFTensor0N();
426 auto t_coords = getFTensor1CoordsAtGaussPts();
427 commonDataPtr->plasticTauJumpPtr->resize(nb_integration_pts, false);
428 commonDataPtr->plasticStrainJumpPtr->resize(6, nb_integration_pts, false);
429 auto t_plastic_strain_jump =
430 getFTensor2SymmetricFromMat<3>(*commonDataPtr->plasticStrainJumpPtr);
431 auto t_tau_jump = getFTensor0FromVec(*commonDataPtr->plasticTauJumpPtr);
432
433 for (size_t gg = 0; gg != nb_integration_pts; ++gg) {
434 double alpha = getMeasure() * t_w;
435
436 // Tensor2_symmetric<PackPtr<double *, 6>, 3> t_nf{&nf[0], &nf[1],
437 // &nf[2],
438 // &nf[3], &nf[4],
439 // &nf[5]};
440 // this is just work in progress
441 t_plastic_strain_jump(i, j) =
442 t_plastic_strain_l(i, j) - t_plastic_strain_r(i, j);
443 t_tau_jump = t_tau_l - t_tau_r;
444
445 // size_t bb = 0;
446 // for (; bb != nb_dofs / 6; ++bb) {
447 // t_nf(i, j) += alpha * t_base;
448 // ++t_base;
449 // ++t_nf;
450 // }
451
452 // for (; bb < nb_base_functions; ++bb)
453 // ++t_base;
454
455 ++t_plastic_strain_jump;
456 ++t_tau_jump;
457 ++t_plastic_strain_l;
458 ++t_plastic_strain_r;
459 ++t_tau_l;
460 ++t_tau_r;
461 ++t_w;
462 }
463
464 // CHKERR VecSetValues(getTSf(), data, nf.data(), ADD_VALUES);
465 }
466
468}
469
471 const std::string field_name, boost::shared_ptr<CommonData> common_data_ptr)
473 commonDataPtr(common_data_ptr) {
474 std::fill(&doEntities[MBVERTEX], &doEntities[MBMAXTYPE], false);
475 doEntities[MBTET] = true;
476 doEntities[MBHEX] = true;
477}
478
480 EntData &data) {
482 const int nb_gauss_pts = getGaussPts().size2();
483 const int nb_base_functions = data.getN().size2();
484 // commonDataPtr->mStressPtr->resize(6, nb_gauss_pts);
485 // auto &t_D = commonDataPtr->tD;
486 // auto t_strain =
487 // getFTensor2SymmetricFromMat<3>(*(commonDataPtr->mStrainPtr)); auto t_stress
488 // = getFTensor2SymmetricFromMat<3>(*(commonDataPtr->mStressPtr));
489 // FIXME: this is wrong, resize only on zero type
490 const size_t nb_dofs = data.getIndices().size();
491 // data.getFieldData().size() == 0
492 // if (side == getFaceSideNumber()) {
493 if (nb_dofs) {
494
495 // for testing //FIXME:
496 MatrixDouble diff = getCoordsAtGaussPts() - getFaceCoordsAtGaussPts();
497 const double eps = 1e-12;
498 if (norm_inf(diff) > eps)
499 SETERRQ(PETSC_COMM_WORLD, MOFEM_ATOM_TEST_INVALID,
500 "coordinates at integration pts are different");
501
502 const size_t nb_gauss_pts = data.getN().size1();
503
504 auto base_function = data.getFTensor0N();
505 // MatrixDouble *ep_elem_data = nullptr;
506 int test_nb_in_loop = getFEMethod()->nInTheLoop;
507 int test_face_sense = getFaceSense();
508 MatrixDouble &mat_ep = commonDataPtr->plasticStrainSideMap[getFaceSense()];
509 mat_ep.resize(6, nb_gauss_pts, false);
510 mat_ep.clear();
511
512 // if (getFEMethod()->nInTheLoop == 0)
513 // ep_elem_data = &(commonDataPtr->plasticStrainLeft);
514 // else
515 // ep_elem_data = &(commonDataPtr->plasticStrainRight);
516
517 // mat_ep.resize(6, nb_gauss_pts, false);
518 // mat_ep.clear();
519 auto values_at_gauss_pts = getFTensor2SymmetricFromMat<3>(mat_ep);
520 for (size_t gg = 0; gg != nb_gauss_pts; ++gg) {
521 auto field_data = data.getFTensor2SymmetricFieldData<3>();
522 size_t bb = 0;
523 for (; bb != nb_dofs / 6; ++bb) {
524 values_at_gauss_pts(i, j) += field_data(i, j) * base_function;
525 ++field_data;
526 ++base_function;
527 }
528 for (; bb != nb_base_functions; ++bb)
529 ++base_function;
530 ++values_at_gauss_pts;
531 }
532 }
534}
535
537 const std::string field_name, boost::shared_ptr<CommonData> common_data_ptr)
539 commonDataPtr(common_data_ptr) {
540 std::fill(&doEntities[MBVERTEX], &doEntities[MBMAXTYPE], false);
541 doEntities[MBTET] = true;
542 doEntities[MBHEX] = true;
543}
544
546 EntData &data) {
548
549 // FIXME: this is wrong, resize only on zero type
550 const size_t nb_dofs = data.getIndices().size();
551 // data.getFieldData().size()
552 if (nb_dofs) {
553 int nb_gauss_pts = data.getN().size1();
554 VectorDouble &vec_tau = commonDataPtr->plasticTauSideMap[getFaceSense()];
555 vec_tau.resize(nb_gauss_pts, false);
556 vec_tau.clear();
557
558 for (int gg = 0; gg != nb_gauss_pts; gg++) {
559 vec_tau(gg) = inner_prod(trans(data.getN(gg)), data.getFieldData());
560 }
561 }
563}
564
566 const std::string field_name, moab::Interface &post_proc_mesh,
567 std::vector<EntityHandle> &map_gauss_pts,
568 boost::shared_ptr<CommonData> common_data_ptr)
570 postProcMesh(post_proc_mesh), mapGaussPts(map_gauss_pts),
571 commonDataPtr(common_data_ptr) {
572 // Operator is only executed for vertices
573 std::fill(&doEntities[MBVERTEX], &doEntities[MBMAXTYPE], false);
574 doEntities[MBTRI] = doEntities[MBQUAD] = true;
575}
576
577//! [Postprocessing]
579 EntData &data) {
581 // if (skinEnts.find(getFEEntityHandle()) != skinEnts.end())
582 // MoFEMFunctionReturnHot(0);
583 const size_t nb_gauss_pts = getGaussPts().size2();
584 std::array<double, 9> def;
585 std::fill(def.begin(), def.end(), 0);
586
587 auto get_tag = [&](const std::string name, size_t size) {
588 Tag th;
589 CHKERR postProcMesh.tag_get_handle(name.c_str(), size, MB_TYPE_DOUBLE, th,
590 MB_TAG_CREAT | MB_TAG_SPARSE,
591 def.data());
592 return th;
593 };
594
595 MatrixDouble3by3 mat(3, 3);
596
597 auto set_matrix_3d = [&](auto &t) -> MatrixDouble3by3 & {
598 mat.clear();
599 for (size_t r = 0; r != 3; ++r)
600 for (size_t c = 0; c != 3; ++c)
601 mat(r, c) = t(r, c);
602 return mat;
603 };
604 auto set_vector = [&](auto &t) -> MatrixDouble3by3 & {
605 mat.clear();
606 for (size_t r = 0; r != 3; ++r)
607 mat(0, r) = t(r);
608 return mat;
609 };
610 auto set_scalar = [&](auto t) -> MatrixDouble3by3 & {
611 mat.clear();
612 mat(0, 0) = t;
613 return mat;
614 };
615
616 auto set_tag = [&](auto th, auto gg, MatrixDouble3by3 &mat) {
617 return postProcMesh.tag_set_data(th, &mapGaussPts[gg], 1,
618 &*mat.data().begin());
619 };
620
621 auto th_rot_vel = get_tag("ROT_VELOCITY", 3);
622 auto th_tau_jump = get_tag("TAU_JUMP", 1);
623 auto th_plastic_jump = get_tag("PLASTIC_STRAIN_JUMP", 9);
624 auto t_coords = getFTensor1CoordsAtGaussPts();
625 auto t_ep_jump =
626 getFTensor2SymmetricFromMat<3>(*(commonDataPtr->plasticStrainJumpPtr));
627 auto t_tau_jump = getFTensor0FromVec(*(commonDataPtr->plasticTauJumpPtr));
628
629 for (size_t gg = 0; gg != nb_gauss_pts; ++gg) {
630
631 Tensor1<double, 3> velocity;
632 velocity(i) = (*cache).Omega(i, j) * t_coords(j);
633
634 CHKERR set_tag(th_tau_jump, gg, set_scalar(t_tau_jump));
635 CHKERR set_tag(th_plastic_jump, gg, set_matrix_3d(t_ep_jump));
636 CHKERR set_tag(th_rot_vel, gg, set_vector(velocity));
637 ++t_coords;
638 ++t_ep_jump;
639 ++t_tau_jump;
640 }
641
643}
644
645} // namespace OpRotatingFrameTools
static const double eps
#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 ...
@ MOFEM_ATOM_TEST_INVALID
Definition definitions.h:40
@ MOFEM_DATA_INCONSISTENCY
Definition definitions.h:31
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
FTensor::Index< 'n', SPACE_DIM > n
FTensor::Index< 'i', SPACE_DIM > i
const double c
speed of light (cm/ns)
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
Tensors class implemented by Walter Landry.
Definition FTensor.hpp:51
Tensor2_Expr< Kronecker_Delta< T >, T, Dim0, Dim1, i, j > kronecker_delta(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
Rank 2.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
constexpr double t
plate stiffness
Definition plate.cpp:59
constexpr auto field_name
chrono::high_resolution_clock::time_point stop
chrono::high_resolution_clock::time_point start
bool sYmm
If true assume that matrix is symmetric structure.
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.
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N(const FieldApproximationBase base)
Get base function as Tensor0.
FTensor::Tensor2_symmetric< FTensor::PackPtr< double *,(Tensor_Dim *(Tensor_Dim+1))/2 >, Tensor_Dim getFTensor2SymmetricFieldData)()
Return symmetric FTensor rank 2, i.e. matrix from filed data coeffinects.
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
const VectorDouble & getFieldData() const
get dofs values
const VectorInt & getIndices() const
Get global indices of dofs on entity.
MoFEMErrorCode loopSideVolumes(const string fe_name, VolumeElementForcesAndSourcesCoreOnSide &fe_method)
auto getFTensor1CoordsAtGaussPts()
Get coordinates at integration points assuming linear geometry.
double getMeasure() const
get measure of element
MatrixDouble & getGaussPts()
matrix of integration (Gauss) points for Volume Element
OpCalculateJumpOnSkeleton(const std::string field_name, boost::shared_ptr< CommonData > common_data_ptr, boost::shared_ptr< DomainSideEle > side_op_fe)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
OpCalculatePlasticConvRotatingFrame(const std::string field_name, boost::shared_ptr< CommonData > common_data_ptr)
OpPostProcPlasticSkeleton(const std::string field_name, moab::Interface &post_proc_mesh, std::vector< EntityHandle > &map_gauss_pts, boost::shared_ptr< CommonData > common_data_ptr)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
[Postprocessing]
boost::shared_ptr< CommonData > commonDataPtr
std::vector< EntityHandle > & mapGaussPts
moab::Interface & postProcMesh
OpRotatingFrameBoundaryLhs(const std::string row_field_name, const std::string col_field_name, boost::shared_ptr< CommonData > common_data_ptr, boost::shared_ptr< DomainSideEle > side_op_fe)
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
Operator for bi-linear form, usually to calculate values on left hand side.
OpRotatingFrameBoundaryRhs(const std::string field_name, boost::shared_ptr< CommonData > common_data_ptr, boost::shared_ptr< DomainSideEle > side_op_fe)
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &data)
MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data)
[Stiffness]
OpRotatingFrameLhs(const std::string row_field_name, const std::string col_field_name, boost::shared_ptr< CommonData > common_data_ptr)
boost::shared_ptr< CommonData > commonDataPtr
OpRotatingFrameRhs(const std::string field_name, boost::shared_ptr< CommonData > common_data_ptr)
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &data)
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data)
OpVolumeSideAssembleRowData(const std::string field_name, boost::shared_ptr< CommonData > common_data_ptr)
OpVolumeSideCalculateEP(const std::string field_name, boost::shared_ptr< CommonData > common_data_ptr)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
OpVolumeSideCalculateTAU(const std::string field_name, boost::shared_ptr< CommonData > common_data_ptr)