v0.13.1
Loading...
Searching...
No Matches
PostProcBrokenMeshInMoabBase.hpp
Go to the documentation of this file.
1/**
2 * @file PostProcBrokenMeshInMoabBase.hpp
3 * @author your name (you@domain.com)
4 * @brief
5 * @version 0.1
6 * @date 2022-09-14
7 *
8 * @copyright Copyright (c) 2022
9 *
10 */
11
12#ifndef __POSTPROCBROKENMESHINMOABBASE_HPP
13#define __POSTPROCBROKENMESHINMOABBASE_HPP
14
15namespace MoFEM {
16
17/**
18 * Each element is subdivided on smaller elements, i.e. a reference mesh on
19 * single element is created. Nodes of such reference mesh are used as
20 * integration points at which field values are calculated and to
21 * each node a "moab" tag is attached to store those values.
22 */
24
25 std::vector<MatrixDouble>
26 levelShapeFunctions; //< values of shape functions on refences element
27 //refinement levels
28 std::vector<MatrixDouble>
29 levelGaussPtsOnRefMesh; //< gauss points at refinement levels
30 std::vector<ublas::matrix<int>> levelRef; //< connectivity at refinement level
31
32 EntityHandle startingVertEleHandle; //< starting handle of first vertex
33 std::vector<double *> verticesOnEleArrays; //< array of vertices coordinate
34 EntityHandle startingEleHandle; //< starting handle of first element
35 EntityHandle *eleConn; //< array of elements connectivity
36
37 int countEle; //< count elements
38 int countVertEle; //< count vertices on the mesh
39
40 int nbVertices; //< numer of vertices on the mesh
41 int nbEles; //< number of elements on the mesh
42
44 virtual ~PostProcGenerateRefMeshBase() = default;
45 MoFEMErrorCode getOptions(); //< get options from command
46
48
49 PetscBool hoNodes; //< if true mid nodes are added
50 int defMaxLevel; //< default max number of refinement levels
51 std::string optPrefix; //< prefix for options
52};
53
55 boost::shared_ptr<PostProcGenerateRefMeshBase>;
56
57/**
58 * @brief Element for postprocessing. Uses MoAB to generate post-processing
59 * mesh.
60 *
61 * @tparam T Finite Element Implementation
62 */
63template <EntityType T> struct PostProcGenerateRefMesh;
64
65template <>
68 MoFEMErrorCode generateReferenceElementMesh();
69};
70
71template <>
74 MoFEMErrorCode generateReferenceElementMesh();
75};
76
77template <>
80 MoFEMErrorCode generateReferenceElementMesh();
81};
82
83template <>
86 MoFEMErrorCode generateReferenceElementMesh();
87};
88
89template <>
92 MoFEMErrorCode generateReferenceElementMesh();
93};
94
95template <typename E> struct PostProcBrokenMeshInMoabBase : public E {
96
99
100 std::vector<EntityHandle> mapGaussPts;
102
103 /**
104 * @brief Get vector of vectors associated to integration points
105 *
106 * @return std::vector<EntityHandle>&
107 */
108 inline auto &getMapGaussPts();
109
110 /**
111 * @brief Get postprocessing mesh
112 *
113 * @return moab::Interface&
114 */
115 inline auto &getPostProcMesh();
116
117 /**
118 * @brief Get postprocessing elements
119 *
120 * @return auto&
121 */
122 inline auto &getPostProcElements();
123
124 /**
125 * \brief wrote results in (MOAB) format, use "file_name.h5m"
126 * @param file_name file name (should always end with .h5m)
127 * @return error code
128
129 * \ingroup mofem_fs_post_proc
130
131 */
132 MoFEMErrorCode writeFile(const std::string file_name);
133
134 /**
135 * @brief Generate vertices and elements
136 *
137 * @return MoFEMErrorCode
138 */
140
142
144
145protected:
146 int getRule(int order);
147 int getMaxLevel() const;
148
149 moab::Core coreMesh;
150 moab::Interface &postProcMesh;
151
152 std::map<EntityType, PostProcGenerateRefMeshPtr> refElementsMap;
153};
154
155template <typename E> int PostProcBrokenMeshInMoabBase<E>::getMaxLevel() const {
156 auto get_element_max_dofs_order = [&]() {
157 int max_order = 0;
158 auto dofs_vec = E::getDataVectorDofsPtr();
159 for (auto &dof : *dofs_vec) {
160 const int dof_order = dof->getDofOrder();
161 max_order = (max_order < dof_order) ? dof_order : max_order;
162 };
163 return max_order;
164 };
165 const auto dof_max_order = get_element_max_dofs_order();
166 return (dof_max_order > 0) ? (dof_max_order - 1) / 2 : 0;
167};
168
169template <typename E>
171 MoFEM::Interface &m_field)
172 : E(m_field), postProcMesh(coreMesh) {}
173
174template <typename E>
176 ParallelComm *pcomm_post_proc_mesh =
177 ParallelComm::get_pcomm(&postProcMesh, MYPCOMM_INDEX);
178 if (pcomm_post_proc_mesh != NULL)
179 delete pcomm_post_proc_mesh;
180}
181
182template <typename E> int PostProcBrokenMeshInMoabBase<E>::getRule(int order) {
183 return -1;
184};
185
186template <typename E>
189
190 auto type = type_from_handle(this->getFEEntityHandle());
191
193
194 try {
195 ref_ele = refElementsMap.at(type);
196 } catch (const out_of_range &e) {
197 SETERRQ1(
198 PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
199 "Generation of reference elements for type <%s> is not implemented",
200 moab::CN::EntityTypeName(type));
201 }
202
203 auto set_gauss_pts = [&](auto &level_gauss_pts_on_ref_mesh, auto &level_ref,
204 auto &level_shape_functions,
205
206 auto start_vert_handle, auto start_ele_handle,
207 auto &verts_array, auto &conn, auto &ver_count,
208 auto &ele_count
209
210 ) {
212
213 size_t level = getMaxLevel();
214 level = std::min(level, level_gauss_pts_on_ref_mesh.size() - 1);
215
216 auto &level_ref_gauss_pts = level_gauss_pts_on_ref_mesh[level];
217 auto &level_ref_ele = level_ref[level];
218 auto &shape_functions = level_shape_functions[level];
219 E::gaussPts.resize(level_ref_gauss_pts.size1(), level_ref_gauss_pts.size2(),
220 false);
221 noalias(E::gaussPts) = level_ref_gauss_pts;
222
223 const auto fe_ent = E::numeredEntFiniteElementPtr->getEnt();
224 auto get_fe_coords = [&]() {
225 const EntityHandle *conn;
226 int num_nodes;
228 E::mField.get_moab().get_connectivity(fe_ent, conn, num_nodes, true),
229 "error get connectivity");
230 VectorDouble coords(num_nodes * 3);
232 E::mField.get_moab().get_coords(conn, num_nodes, &*coords.begin()),
233 "error get coordinates");
234 return coords;
235 };
236
237 auto coords = get_fe_coords();
238
239 const int num_nodes = level_ref_gauss_pts.size2();
240 mapGaussPts.resize(level_ref_gauss_pts.size2());
241
244 &*shape_functions.data().begin());
246 &verts_array[0][ver_count], &verts_array[1][ver_count],
247 &verts_array[2][ver_count]);
248 for (int gg = 0; gg != num_nodes; ++gg, ++ver_count) {
249
250 mapGaussPts[gg] = start_vert_handle + ver_count;
251
252 auto set_float_precision = [](const double x) {
253 if (std::abs(x) < std::numeric_limits<float>::epsilon())
254 return 0.;
255 else
256 return x;
257 };
258
259 t_coords(i) = 0;
260 auto t_ele_coords = getFTensor1FromArray<3, 3>(coords);
261 for (int nn = 0; nn != CN::VerticesPerEntity(type); ++nn) {
262 t_coords(i) += t_n * t_ele_coords(i);
263 ++t_ele_coords;
264 ++t_n;
265 }
266
267 for (auto ii : {0, 1, 2})
268 t_coords(ii) = set_float_precision(t_coords(ii));
269
270 ++t_coords;
271 }
272
273 Tag th;
274 int def_in_the_loop = -1;
275 CHKERR postProcMesh.tag_get_handle("NB_IN_THE_LOOP", 1, MB_TYPE_INTEGER, th,
276 MB_TAG_CREAT | MB_TAG_SPARSE,
277 &def_in_the_loop);
278
279 postProcElements.clear();
280 const int num_el = level_ref_ele.size1();
281 const int num_nodes_on_ele = level_ref_ele.size2();
282 auto start_e = start_ele_handle + ele_count;
283 postProcElements = Range(start_e, start_e + num_el - 1);
284 for (auto tt = 0; tt != level_ref_ele.size1(); ++tt, ++ele_count) {
285 for (int nn = 0; nn != num_nodes_on_ele; ++nn) {
286 conn[num_nodes_on_ele * ele_count + nn] =
287 mapGaussPts[level_ref_ele(tt, nn)];
288 }
289 }
290
291 const int n_in_the_loop = E::nInTheLoop;
292 CHKERR postProcMesh.tag_clear_data(th, postProcElements, &n_in_the_loop);
293
295 };
296
297 CHKERR set_gauss_pts(
298
299 ref_ele->levelGaussPtsOnRefMesh, ref_ele->levelRef,
300 ref_ele->levelShapeFunctions,
301
302 ref_ele->startingVertEleHandle, ref_ele->startingEleHandle,
303 ref_ele->verticesOnEleArrays, ref_ele->eleConn, ref_ele->countVertEle,
304 ref_ele->countEle
305
306 );
307
309};
310
311template <typename E>
313 moab::Interface &moab = coreMesh;
315
316 ParallelComm *pcomm_post_proc_mesh =
317 ParallelComm::get_pcomm(&moab, MYPCOMM_INDEX);
318 if (pcomm_post_proc_mesh != NULL)
319 delete pcomm_post_proc_mesh;
320
321 CHKERR postProcMesh.delete_mesh();
322
323 auto get_ref_ele = [&](const EntityType type) {
324 PostProcGenerateRefMeshPtr ref_ele_ptr;
325
326 try {
327
328 ref_ele_ptr = refElementsMap.at(type);
329
330 } catch (const out_of_range &e) {
331
332 switch (type) {
333 case MBTET:
334 ref_ele_ptr = boost::make_shared<PostProcGenerateRefMesh<MBTET>>();
335 break;
336 case MBHEX:
337 ref_ele_ptr = boost::make_shared<PostProcGenerateRefMesh<MBHEX>>();
338 break;
339 case MBTRI:
340 ref_ele_ptr = boost::make_shared<PostProcGenerateRefMesh<MBTRI>>();
341 break;
342 case MBQUAD:
343 ref_ele_ptr = boost::make_shared<PostProcGenerateRefMesh<MBQUAD>>();
344 break;
345 case MBEDGE:
346 ref_ele_ptr = boost::make_shared<PostProcGenerateRefMesh<MBEDGE>>();
347 break;
348 default:
349 MOFEM_LOG("SELF", Sev::error)
350 << "Generation of reference elements for type < "
351 << moab::CN::EntityTypeName(type) << " > is not implemented";
352 CHK_THROW_MESSAGE(MOFEM_NOT_IMPLEMENTED, "Element not implemented");
353 }
354
355 CHK_THROW_MESSAGE(ref_ele_ptr->generateReferenceElementMesh(),
356 "Error when generating reference element");
357 }
358
359 refElementsMap[type] = ref_ele_ptr;
360
361 return ref_ele_ptr;
362 };
363
364 auto fe_ptr = this->problemPtr->numeredFiniteElementsPtr;
365
366 auto miit =
367 fe_ptr->template get<Composite_Name_And_Part_mi_tag>().lower_bound(
368 boost::make_tuple(this->getFEName(), this->getLoFERank()));
369 auto hi_miit =
370 fe_ptr->template get<Composite_Name_And_Part_mi_tag>().upper_bound(
371 boost::make_tuple(this->getFEName(), this->getHiFERank()));
372
373 const int number_of_ents_in_the_loop = this->getLoopSize();
374 if (std::distance(miit, hi_miit) != number_of_ents_in_the_loop) {
375 SETERRQ(E::mField.get_comm(), MOFEM_DATA_INCONSISTENCY,
376 "Wrong size of indicies. Inconsistent size number of iterated "
377 "elements iterated by problem and from range.");
378 }
379
380 for (auto &m : refElementsMap) {
381 m.second->nbVertices = 0;
382 m.second->nbEles = 0;
383 m.second->countEle = 0;
384 m.second->countVertEle = 0;
385 }
386
387 for (; miit != hi_miit; ++miit) {
388 auto type = (*miit)->getEntType();
389 auto ref_ele = get_ref_ele(type);
390
391 // Set pointer to element. So that getDataVectorDofsPtr in getMaxLevel
392 // can work
393 E::numeredEntFiniteElementPtr = *miit;
394 bool add = true;
395 if (E::exeTestHook) {
396 add = E::exeTestHook(this);
397 }
398
399 if (add) {
400 size_t level = getMaxLevel();
401 level = std::min(level, ref_ele->levelGaussPtsOnRefMesh.size() - 1);
402 ref_ele->nbVertices += ref_ele->levelGaussPtsOnRefMesh[level].size2();
403 ref_ele->nbEles += ref_ele->levelRef[level].size1();
404 }
405 }
406
407 auto alloc_vertices_and_elements_on_post_proc_mesh = [&]() {
409
410 ReadUtilIface *iface;
411 CHKERR postProcMesh.query_interface(iface);
412
413 for (auto &m : refElementsMap) {
414 if (m.second->nbEles) {
415 CHKERR iface->get_node_coords(3, m.second->nbVertices, 0,
416 m.second->startingVertEleHandle,
417 m.second->verticesOnEleArrays);
418 CHKERR iface->get_element_connect(
419 m.second->nbEles, m.second->levelRef[0].size2(), m.first, 0,
420 m.second->startingEleHandle, m.second->eleConn);
421
422 m.second->countEle = 0;
423 m.second->countVertEle = 0;
424 }
425 }
426
428 };
429
430 CHKERR alloc_vertices_and_elements_on_post_proc_mesh();
431
433}
434
435template <typename E>
438
439 auto update_elements = [&]() {
440 ReadUtilIface *iface;
441 CHKERR this->postProcMesh.query_interface(iface);
443
444 for (auto &m : refElementsMap) {
445 if (m.second->nbEles) {
446 MOFEM_TAG_AND_LOG("SELF", Sev::noisy, "PostProc")
447 << "Update < " << moab::CN::EntityTypeName(m.first)
448 << m.second->countEle;
449 CHKERR iface->update_adjacencies(
450 m.second->startingEleHandle, m.second->countEle,
451 m.second->levelRef[0].size2(), m.second->eleConn);
452 }
453 }
454
456 };
457
458 auto resolve_shared_ents = [&]() {
460 auto get_lower_dimension = [&]() {
461 int min_dim = 3;
462 for (auto &m : refElementsMap) {
463 const int dim = moab::CN::Dimension(m.first);
464 if (m.second->nbEles) {
465 min_dim = std::min(dim, min_dim);
466 }
467 }
468 return min_dim;
469 };
470
471 auto remove_obsolete_entities = [&](auto &&min_dim) {
472 Range edges;
473 CHKERR postProcMesh.get_entities_by_type(0, MBEDGE, edges, false);
474
475 Range faces;
476 CHKERR postProcMesh.get_entities_by_dimension(0, 2, faces, false);
477
478 Range vols;
479 CHKERR postProcMesh.get_entities_by_dimension(0, 3, vols, false);
480
481 Range ents;
482
483 if (min_dim > 1)
484 CHKERR postProcMesh.delete_entities(edges);
485 else
486 ents.merge(edges);
487
488 if (min_dim > 2)
489 CHKERR postProcMesh.delete_entities(faces);
490 else
491 ents.merge(faces);
492
493 ents.merge(vols);
494
495 return ents;
496 };
497
498 auto ents = remove_obsolete_entities(get_lower_dimension());
499
500 ParallelComm *pcomm_post_proc_mesh =
501 ParallelComm::get_pcomm(&(postProcMesh), MYPCOMM_INDEX);
502 if (pcomm_post_proc_mesh == NULL) {
503 // wrapRefMeshComm =
504 // boost::make_shared<WrapMPIComm>(T::mField.get_comm(), false);
505 pcomm_post_proc_mesh = new ParallelComm(
506 &(postProcMesh),
507 PETSC_COMM_WORLD /*(T::wrapRefMeshComm)->get_comm()*/);
508 }
509
510 int rank = E::mField.get_comm_rank();
511 CHKERR postProcMesh.tag_clear_data(pcomm_post_proc_mesh->part_tag(), ents,
512 &rank);
513
514 CHKERR pcomm_post_proc_mesh->resolve_shared_ents(0);
515
517 };
518
519 CHKERR update_elements();
520 CHKERR resolve_shared_ents();
521
523}
524
526 return mapGaussPts;
527}
528
530 return postProcMesh;
531}
532
533template <typename E>
535 return postProcElements;
536}
537
538template <typename E>
540PostProcBrokenMeshInMoabBase<E>::writeFile(const std::string file_name) {
542 ParallelComm *pcomm_post_proc_mesh =
543 ParallelComm::get_pcomm(&postProcMesh, MYPCOMM_INDEX);
544 if (pcomm_post_proc_mesh == NULL)
545 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
546 "ParallelComm not allocated");
547 CHKERR postProcMesh.write_file(file_name.c_str(), "MOAB",
548 "PARALLEL=WRITE_PART");
550};
551
552template <typename E> struct PostProcBrokenMeshInMoab;
553
554template <>
556 : public PostProcBrokenMeshInMoabBase<VolumeElementForcesAndSourcesCore> {
558 VolumeElementForcesAndSourcesCore>::PostProcBrokenMeshInMoabBase;
559};
560
561template <>
563 : public PostProcBrokenMeshInMoabBase<FaceElementForcesAndSourcesCore> {
565 FaceElementForcesAndSourcesCore>::PostProcBrokenMeshInMoabBase;
566};
567
568template <>
570 : public PostProcBrokenMeshInMoabBase<EdgeElementForcesAndSourcesCore> {
572 EdgeElementForcesAndSourcesCore>::PostProcBrokenMeshInMoabBase;
573};
574
575/**
576 * @brief Post post-proc data at points from hash maps
577 *
578 * @tparam DIM1 dimension of vector in data_map_vec and first column of
579 * data_map_may
580 * @tparam DIM2 dimension of second column in data_map_mat
581 */
582template <int DIM1, int DIM2>
584
585 using DataMapVec = std::map<std::string, boost::shared_ptr<VectorDouble>>;
586 using DataMapMat = std::map<std::string, boost::shared_ptr<MatrixDouble>>;
587
588 /**
589 * @brief Construct a new OpPostProcMapInMoab object
590 *
591 * @param post_proc_mesh postprocessing mesh
592 * @param map_gauss_pts map of gauss points to nodes of postprocessing mesh
593 * @param data_map_scalar hash map of scalar values (string is name of the
594 * tag)
595 * @param data_map_vec hash map of vector values
596 * @param data_map_mat hash map of second order tensor values
597 * @param data_symm_map_mat hash map of symmetric second order tensor values
598 */
599 OpPostProcMapInMoab(moab::Interface &post_proc_mesh,
600 std::vector<EntityHandle> &map_gauss_pts,
601 DataMapVec data_map_scalar, DataMapMat data_map_vec,
602 DataMapMat data_map_mat, DataMapMat data_symm_map_mat)
605 postProcMesh(post_proc_mesh), mapGaussPts(map_gauss_pts),
606 dataMapScalar(data_map_scalar), dataMapVec(data_map_vec),
607 dataMapMat(data_map_mat), dataMapSymmMat(data_symm_map_mat) {
608 // Operator is only executed for vertices
609 std::fill(&doEntities[MBEDGE], &doEntities[MBMAXTYPE], false);
610 }
611 MoFEMErrorCode doWork(int side, EntityType type,
613
614private:
615 moab::Interface &postProcMesh;
616 std::vector<EntityHandle> &mapGaussPts;
621};
622
623template <int DIM1, int DIM2>
628
629 std::array<double, 9> def;
630 std::fill(def.begin(), def.end(), 0);
631
632 auto get_tag = [&](const std::string name, size_t size) {
633 Tag th;
634 CHKERR postProcMesh.tag_get_handle(name.c_str(), size, MB_TYPE_DOUBLE, th,
635 MB_TAG_CREAT | MB_TAG_SPARSE,
636 def.data());
637 return th;
638 };
639
640 MatrixDouble3by3 mat(3, 3);
641
642 auto set_vector_3d = [&](auto &t) -> MatrixDouble3by3 & {
643 mat.clear();
644 for (size_t r = 0; r != DIM1; ++r)
645 mat(0, r) = t(r);
646 return mat;
647 };
648
649 auto set_matrix_3d = [&](auto &t) -> MatrixDouble3by3 & {
650 mat.clear();
651 for (size_t r = 0; r != DIM1; ++r)
652 for (size_t c = 0; c != DIM2; ++c)
653 mat(r, c) = t(r, c);
654 return mat;
655 };
656
657 auto set_matrix_symm_3d = [&](auto &t) -> MatrixDouble3by3 & {
658 mat.clear();
659 for (size_t r = 0; r != DIM1; ++r)
660 for (size_t c = 0; c != DIM1; ++c)
661 mat(r, c) = t(r, c);
662 return mat;
663 };
664
665 auto set_scalar = [&](auto t) -> MatrixDouble3by3 & {
666 mat.clear();
667 mat(0, 0) = t;
668 return mat;
669 };
670
671 auto set_float_precision = [](const double x) {
672 if (std::abs(x) < std::numeric_limits<float>::epsilon())
673 return 0.;
674 else
675 return x;
676 };
677
678 auto set_tag = [&](auto th, auto gg, MatrixDouble3by3 &mat) {
679 for (auto &v : mat.data())
680 v = set_float_precision(v);
681 return postProcMesh.tag_set_data(th, &mapGaussPts[gg], 1,
682 &*mat.data().begin());
683 };
684
685 for (auto &m : dataMapScalar) {
686 auto th = get_tag(m.first, 1);
687 auto t_scl = getFTensor0FromVec(*m.second);
688 auto nb_integration_pts = getGaussPts().size2();
689 size_t gg = 0;
690 for (int gg = 0; gg != nb_integration_pts; ++gg) {
691 CHKERR set_tag(th, gg, set_scalar(t_scl));
692 ++t_scl;
693 }
694 }
695
696 for (auto &m : dataMapVec) {
697 auto th = get_tag(m.first, 3);
698 auto t_vec = getFTensor1FromMat<DIM1>(*m.second);
699 auto nb_integration_pts = getGaussPts().size2();
700 size_t gg = 0;
701 for (int gg = 0; gg != nb_integration_pts; ++gg) {
702 CHKERR set_tag(th, gg, set_vector_3d(t_vec));
703 ++t_vec;
704 }
705 }
706
707 for (auto &m : dataMapMat) {
708 auto th = get_tag(m.first, 9);
709 auto t_mat = getFTensor2FromMat<DIM1, DIM2>(*m.second);
710 auto nb_integration_pts = getGaussPts().size2();
711 size_t gg = 0;
712 for (int gg = 0; gg != nb_integration_pts; ++gg) {
713 CHKERR set_tag(th, gg, set_matrix_3d(t_mat));
714 ++t_mat;
715 }
716 }
717
718 for (auto &m : dataMapSymmMat) {
719 auto th = get_tag(m.first, 9);
720 auto t_mat = getFTensor2SymmetricFromMat<DIM1>(*m.second);
721 auto nb_integration_pts = getGaussPts().size2();
722 size_t gg = 0;
723 for (int gg = 0; gg != nb_integration_pts; ++gg) {
724 CHKERR set_tag(th, gg, set_matrix_symm_3d(t_mat));
725 ++t_mat;
726 }
727 }
728
730}
731} // namespace MoFEM
732
733#endif //__POSTPROCBROKENMESHINMOABBASE_HPP
ForcesAndSourcesCore::UserDataOperator UserDataOperator
#define MOFEM_TAG_AND_LOG(channel, severity, tag)
Tag and log in channel.
Definition: LogManager.hpp:345
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
Definition: definitions.h:595
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
@ NOSPACE
Definition: definitions.h:83
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:215
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define CHK_MOAB_THROW(err, msg)
Check error code of MoAB function and throw MoFEM exception.
Definition: definitions.h:576
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
#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
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
const int dim
FTensor::Index< 'm', SPACE_DIM > m
MoFEMErrorCode writeFile(const std::string file_name)
wrote results in (MOAB) format, use "file_name.h5m"
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:301
FTensor::Index< 'i', SPACE_DIM > i
const double c
speed of light (cm/ns)
const double v
phase velocity of light in medium (cm/ns)
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
auto type_from_handle(const EntityHandle h)
get type from entity handle
Definition: Templates.hpp:1592
boost::shared_ptr< PostProcGenerateRefMeshBase > PostProcGenerateRefMeshPtr
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
Definition: Templates.hpp:135
constexpr double t
plate stiffness
Definition: plate.cpp:59
std::array< bool, MBMAXTYPE > doEntities
If true operator is executed for entity.
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
@ OPSPACE
operator do Work is execute on space data
structure to get information form mofem into EntitiesFieldData
Post post-proc data at points from hash maps.
DataMapMat dataMapMat
DataMapMat dataMapVec
DataMapMat dataMapSymmMat
std::map< std::string, boost::shared_ptr< VectorDouble > > DataMapVec
std::map< std::string, boost::shared_ptr< MatrixDouble > > DataMapMat
moab::Interface & postProcMesh
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
std::vector< EntityHandle > & mapGaussPts
DataMapVec dataMapScalar
OpPostProcMapInMoab(moab::Interface &post_proc_mesh, std::vector< EntityHandle > &map_gauss_pts, DataMapVec data_map_scalar, DataMapMat data_map_vec, DataMapMat data_map_mat, DataMapMat data_symm_map_mat)
Construct a new OpPostProcMapInMoab object.
auto & getPostProcElements()
Get postprocessing elements.
MoFEMErrorCode preProcess()
Generate vertices and elements.
moab::Core coreMesh
MoFEMErrorCode setGaussPts(int order)
moab::Interface & postProcMesh
PostProcBrokenMeshInMoabBase(MoFEM::Interface &m_field)
auto & getMapGaussPts()
Get vector of vectors associated to integration points.
auto & getPostProcMesh()
Get postprocessing mesh.
MoFEMErrorCode postProcess()
std::vector< EntityHandle > mapGaussPts
int getMaxLevel() const
Range postProcElements
virtual ~PostProcBrokenMeshInMoabBase()
std::map< EntityType, PostProcGenerateRefMeshPtr > refElementsMap
int getRule(int order)
EntityHandle startingEleHandle
int countVertEle
virtual MoFEMErrorCode generateReferenceElementMesh()=0
virtual ~PostProcGenerateRefMeshBase()=default
EntityHandle * eleConn
std::vector< ublas::matrix< int > > levelRef
std::vector< double * > verticesOnEleArrays
PetscBool hoNodes
int defMaxLevel
std::vector< MatrixDouble > levelShapeFunctions
PostProcGenerateRefMeshBase()
int countEle
int nbVertices
std::vector< MatrixDouble > levelGaussPtsOnRefMesh
int nbEles
MoFEMErrorCode getOptions()
std::string optPrefix
EntityHandle startingVertEleHandle
Element for postprocessing. Uses MoAB to generate post-processing mesh.