v0.13.2
Loading...
Searching...
No Matches
Public Member Functions | Public Attributes | Protected Member Functions | Protected Attributes | List of all members
MoFEM::PostProcBrokenMeshInMoabBase< E > Struct Template Reference

#include <src/post_proc/PostProcBrokenMeshInMoabBase.hpp>

Inheritance diagram for MoFEM::PostProcBrokenMeshInMoabBase< E >:
[legend]
Collaboration diagram for MoFEM::PostProcBrokenMeshInMoabBase< E >:
[legend]

Public Member Functions

 PostProcBrokenMeshInMoabBase (MoFEM::Interface &m_field)
 
virtual ~PostProcBrokenMeshInMoabBase ()
 
auto & getMapGaussPts ()
 Get vector of vectors associated to integration points. More...
 
auto & getPostProcMesh ()
 Get postprocessing mesh. More...
 
auto & getPostProcElements ()
 Get postprocessing elements. More...
 
MoFEMErrorCode writeFile (const std::string file_name)
 wrote results in (MOAB) format, use "file_name.h5m" More...
 
MoFEMErrorCode preProcess ()
 Generate vertices and elements. More...
 
MoFEMErrorCode postProcess ()
 
MoFEMErrorCode setGaussPts (int order)
 
MoFEMErrorCode setTagsToTransfer (std::vector< Tag > tags_to_transfer)
 

Public Attributes

std::vector< EntityHandlemapGaussPts
 
Range postProcElements
 

Protected Member Functions

int getRule (int order)
 
virtual int getMaxLevel () const
 Determine refinment level based on fields approx ordre. More...
 
virtual MoFEMErrorCode transferTags ()
 

Protected Attributes

moab::Core coreMesh
 
moab::Interface & postProcMesh
 
std::map< EntityType, PostProcGenerateRefMeshPtrrefElementsMap
 
std::vector< Tag > tagsToTransfer
 

Detailed Description

template<typename E>
struct MoFEM::PostProcBrokenMeshInMoabBase< E >

Definition at line 95 of file PostProcBrokenMeshInMoabBase.hpp.

Constructor & Destructor Documentation

◆ PostProcBrokenMeshInMoabBase()

template<typename E >
MoFEM::PostProcBrokenMeshInMoabBase< E >::PostProcBrokenMeshInMoabBase ( MoFEM::Interface m_field)

◆ ~PostProcBrokenMeshInMoabBase()

Definition at line 257 of file PostProcBrokenMeshInMoabBase.hpp.

257 {
258 ParallelComm *pcomm_post_proc_mesh =
259 ParallelComm::get_pcomm(&postProcMesh, MYPCOMM_INDEX);
260 if (pcomm_post_proc_mesh != NULL)
261 delete pcomm_post_proc_mesh;
262}
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:215

Member Function Documentation

◆ getMapGaussPts()

template<typename E >
auto & MoFEM::PostProcBrokenMeshInMoabBase< E >::getMapGaussPts
inline

Get vector of vectors associated to integration points.

Returns
std::vector<EntityHandle>&

Definition at line 608 of file PostProcBrokenMeshInMoabBase.hpp.

608 {
609 return mapGaussPts;
610}
std::vector< EntityHandle > mapGaussPts

◆ getMaxLevel()

template<typename E >
int MoFEM::PostProcBrokenMeshInMoabBase< E >::getMaxLevel
protectedvirtual

Determine refinment level based on fields approx ordre.

level = (order - 1) / 2

Returns
int

Definition at line 179 of file PostProcBrokenMeshInMoabBase.hpp.

179 {
180 auto get_element_max_dofs_order = [&]() {
181 int max_order = 0;
182 auto dofs_vec = E::getDataVectorDofsPtr();
183 for (auto &dof : *dofs_vec) {
184 const int dof_order = dof->getDofOrder();
185 max_order = (max_order < dof_order) ? dof_order : max_order;
186 };
187 return max_order;
188 };
189 const auto dof_max_order = get_element_max_dofs_order();
190 return (dof_max_order > 0) ? (dof_max_order - 1) / 2 : 0;
191};

◆ getPostProcElements()

template<typename E >
auto & MoFEM::PostProcBrokenMeshInMoabBase< E >::getPostProcElements
inline

Get postprocessing elements.

Returns
auto&

Definition at line 617 of file PostProcBrokenMeshInMoabBase.hpp.

◆ getPostProcMesh()

template<typename E >
auto & MoFEM::PostProcBrokenMeshInMoabBase< E >::getPostProcMesh
inline

Get postprocessing mesh.

Returns
moab::Interface&

Definition at line 612 of file PostProcBrokenMeshInMoabBase.hpp.

612 {
613 return postProcMesh;
614}

◆ getRule()

template<typename E >
int MoFEM::PostProcBrokenMeshInMoabBase< E >::getRule ( int  order)
protected

Definition at line 264 of file PostProcBrokenMeshInMoabBase.hpp.

264 {
265 return -1;
266};

◆ postProcess()

template<typename E >
MoFEMErrorCode MoFEM::PostProcBrokenMeshInMoabBase< E >::postProcess

Definition at line 519 of file PostProcBrokenMeshInMoabBase.hpp.

519 {
521
522 auto update_elements = [&]() {
523 ReadUtilIface *iface;
524 CHKERR this->postProcMesh.query_interface(iface);
526
527 for (auto &m : refElementsMap) {
528 if (m.second->nbEles) {
529 MOFEM_TAG_AND_LOG("SELF", Sev::noisy, "PostProc")
530 << "Update < " << moab::CN::EntityTypeName(m.first)
531 << " number of processsed " << m.second->countEle;
532 CHKERR iface->update_adjacencies(
533 m.second->startingEleHandle, m.second->countEle,
534 m.second->levelRef[0].size2(), m.second->eleConn);
535 }
536 }
537
539 };
540
541 auto resolve_shared_ents = [&]() {
543 auto get_lower_dimension = [&]() {
544 int min_dim = 3;
545 for (auto &m : refElementsMap) {
546 const int dim = moab::CN::Dimension(m.first);
547 if (m.second->nbEles) {
548 min_dim = std::min(dim, min_dim);
549 }
550 }
551 return min_dim;
552 };
553
554 auto remove_obsolete_entities = [&](auto &&min_dim) {
555 Range edges;
556 CHKERR postProcMesh.get_entities_by_type(0, MBEDGE, edges, false);
557
558 Range faces;
559 CHKERR postProcMesh.get_entities_by_dimension(0, 2, faces, false);
560
561 Range vols;
562 CHKERR postProcMesh.get_entities_by_dimension(0, 3, vols, false);
563
564 Range ents;
565
566 if (min_dim > 1)
567 CHKERR postProcMesh.delete_entities(edges);
568 else
569 ents.merge(edges);
570
571 if (min_dim > 2)
572 CHKERR postProcMesh.delete_entities(faces);
573 else
574 ents.merge(faces);
575
576 ents.merge(vols);
577
578 return ents;
579 };
580
581 auto ents = remove_obsolete_entities(get_lower_dimension());
582
583 ParallelComm *pcomm_post_proc_mesh =
584 ParallelComm::get_pcomm(&(postProcMesh), MYPCOMM_INDEX);
585 if (pcomm_post_proc_mesh == NULL) {
586 // wrapRefMeshComm =
587 // boost::make_shared<WrapMPIComm>(T::mField.get_comm(), false);
588 pcomm_post_proc_mesh = new ParallelComm(
589 &(postProcMesh),
590 PETSC_COMM_WORLD /*(T::wrapRefMeshComm)->get_comm()*/);
591 }
592
593 int rank = E::mField.get_comm_rank();
594 CHKERR postProcMesh.tag_clear_data(pcomm_post_proc_mesh->part_tag(), ents,
595 &rank);
596
597 CHKERR pcomm_post_proc_mesh->resolve_shared_ents(0);
598
600 };
601
602 CHKERR update_elements();
603 CHKERR resolve_shared_ents();
604
606}
#define MOFEM_TAG_AND_LOG(channel, severity, tag)
Tag and log in channel.
Definition: LogManager.hpp:345
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#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
#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
std::map< EntityType, PostProcGenerateRefMeshPtr > refElementsMap

◆ preProcess()

template<typename E >
MoFEMErrorCode MoFEM::PostProcBrokenMeshInMoabBase< E >::preProcess

Generate vertices and elements.

Returns
MoFEMErrorCode

Definition at line 395 of file PostProcBrokenMeshInMoabBase.hpp.

395 {
396 moab::Interface &moab = coreMesh;
398
399 ParallelComm *pcomm_post_proc_mesh =
400 ParallelComm::get_pcomm(&moab, MYPCOMM_INDEX);
401 if (pcomm_post_proc_mesh != NULL)
402 delete pcomm_post_proc_mesh;
403
404 CHKERR postProcMesh.delete_mesh();
405
406 auto get_ref_ele = [&](const EntityType type) {
407 PostProcGenerateRefMeshPtr ref_ele_ptr;
408
409 try {
410
411 ref_ele_ptr = refElementsMap.at(type);
412
413 } catch (const out_of_range &e) {
414
415 switch (type) {
416 case MBTET:
417 ref_ele_ptr = boost::make_shared<PostProcGenerateRefMesh<MBTET>>();
418 break;
419 case MBHEX:
420 ref_ele_ptr = boost::make_shared<PostProcGenerateRefMesh<MBHEX>>();
421 break;
422 case MBTRI:
423 ref_ele_ptr = boost::make_shared<PostProcGenerateRefMesh<MBTRI>>();
424 break;
425 case MBQUAD:
426 ref_ele_ptr = boost::make_shared<PostProcGenerateRefMesh<MBQUAD>>();
427 break;
428 case MBEDGE:
429 ref_ele_ptr = boost::make_shared<PostProcGenerateRefMesh<MBEDGE>>();
430 break;
431 default:
432 MOFEM_LOG("SELF", Sev::error)
433 << "Generation of reference elements for type < "
434 << moab::CN::EntityTypeName(type) << " > is not implemented";
435 CHK_THROW_MESSAGE(MOFEM_NOT_IMPLEMENTED, "Element not implemented");
436 }
437
438 CHK_THROW_MESSAGE(ref_ele_ptr->generateReferenceElementMesh(),
439 "Error when generating reference element");
440 }
441
442 refElementsMap[type] = ref_ele_ptr;
443
444 return ref_ele_ptr;
445 };
446
447 auto fe_ptr = this->problemPtr->numeredFiniteElementsPtr;
448
449 auto miit =
450 fe_ptr->template get<Composite_Name_And_Part_mi_tag>().lower_bound(
451 boost::make_tuple(this->getFEName(), this->getLoFERank()));
452 auto hi_miit =
453 fe_ptr->template get<Composite_Name_And_Part_mi_tag>().upper_bound(
454 boost::make_tuple(this->getFEName(), this->getHiFERank()));
455
456 const int number_of_ents_in_the_loop = this->getLoopSize();
457 if (std::distance(miit, hi_miit) != number_of_ents_in_the_loop) {
458 SETERRQ(E::mField.get_comm(), MOFEM_DATA_INCONSISTENCY,
459 "Wrong size of indicies. Inconsistent size number of iterated "
460 "elements iterated by problem and from range.");
461 }
462
463 for (auto &m : refElementsMap) {
464 m.second->nbVertices = 0;
465 m.second->nbEles = 0;
466 m.second->countEle = 0;
467 m.second->countVertEle = 0;
468 }
469
470 for (; miit != hi_miit; ++miit) {
471 auto type = (*miit)->getEntType();
472 auto ref_ele = get_ref_ele(type);
473
474 // Set pointer to element. So that getDataVectorDofsPtr in getMaxLevel
475 // can work
476 E::numeredEntFiniteElementPtr = *miit;
477 bool add = true;
478 if (E::exeTestHook) {
479 add = E::exeTestHook(this);
480 }
481
482 if (add) {
483 size_t level = getMaxLevel();
484 level = std::min(level, ref_ele->levelGaussPtsOnRefMesh.size() - 1);
485 ref_ele->nbVertices += ref_ele->levelGaussPtsOnRefMesh[level].size2();
486 ref_ele->nbEles += ref_ele->levelRef[level].size1();
487 }
488 }
489
490 auto alloc_vertices_and_elements_on_post_proc_mesh = [&]() {
492
493 ReadUtilIface *iface;
494 CHKERR postProcMesh.query_interface(iface);
495
496 for (auto &m : refElementsMap) {
497 if (m.second->nbEles) {
498 CHKERR iface->get_node_coords(3, m.second->nbVertices, 0,
499 m.second->startingVertEleHandle,
500 m.second->verticesOnEleArrays);
501 CHKERR iface->get_element_connect(
502 m.second->nbEles, m.second->levelRef[0].size2(), m.first, 0,
503 m.second->startingEleHandle, m.second->eleConn);
504
505 m.second->countEle = 0;
506 m.second->countVertEle = 0;
507 }
508 }
509
511 };
512
513 CHKERR alloc_vertices_and_elements_on_post_proc_mesh();
514
516}
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
Definition: definitions.h:595
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:301
boost::shared_ptr< PostProcGenerateRefMeshBase > PostProcGenerateRefMeshPtr
virtual int getMaxLevel() const
Determine refinment level based on fields approx ordre.

◆ setGaussPts()

template<typename E >
MoFEMErrorCode MoFEM::PostProcBrokenMeshInMoabBase< E >::setGaussPts ( int  order)

Definition at line 269 of file PostProcBrokenMeshInMoabBase.hpp.

269 {
271
272 auto type = type_from_handle(this->getFEEntityHandle());
273
275
276 try {
277 ref_ele = refElementsMap.at(type);
278 } catch (const out_of_range &e) {
279 SETERRQ1(
280 PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
281 "Generation of reference elements for type <%s> is not implemented",
282 moab::CN::EntityTypeName(type));
283 }
284
285 auto set_gauss_pts = [&](auto &level_gauss_pts_on_ref_mesh, auto &level_ref,
286 auto &level_shape_functions,
287
288 auto start_vert_handle, auto start_ele_handle,
289 auto &verts_array, auto &conn, auto &ver_count,
290 auto &ele_count
291
292 ) {
294
295 size_t level = getMaxLevel();
296 level = std::min(level, level_gauss_pts_on_ref_mesh.size() - 1);
297
298 auto &level_ref_gauss_pts = level_gauss_pts_on_ref_mesh[level];
299 auto &level_ref_ele = level_ref[level];
300 auto &shape_functions = level_shape_functions[level];
301 E::gaussPts.resize(level_ref_gauss_pts.size1(), level_ref_gauss_pts.size2(),
302 false);
303 noalias(E::gaussPts) = level_ref_gauss_pts;
304
305 const auto fe_ent = E::numeredEntFiniteElementPtr->getEnt();
306 auto get_fe_coords = [&]() {
307 const EntityHandle *conn;
308 int num_nodes;
310 E::mField.get_moab().get_connectivity(fe_ent, conn, num_nodes, true),
311 "error get connectivity");
312 VectorDouble coords(num_nodes * 3);
314 E::mField.get_moab().get_coords(conn, num_nodes, &*coords.begin()),
315 "error get coordinates");
316 return coords;
317 };
318
319 auto coords = get_fe_coords();
320
321 const int num_nodes = level_ref_gauss_pts.size2();
322 mapGaussPts.resize(level_ref_gauss_pts.size2());
323
326 &*shape_functions.data().begin());
328 &verts_array[0][ver_count], &verts_array[1][ver_count],
329 &verts_array[2][ver_count]);
330 for (int gg = 0; gg != num_nodes; ++gg, ++ver_count) {
331
332 mapGaussPts[gg] = start_vert_handle + ver_count;
333
334 auto set_float_precision = [](const double x) {
335 if (std::abs(x) < std::numeric_limits<float>::epsilon())
336 return 0.;
337 else
338 return x;
339 };
340
341 t_coords(i) = 0;
342 auto t_ele_coords = getFTensor1FromArray<3, 3>(coords);
343 for (int nn = 0; nn != CN::VerticesPerEntity(type); ++nn) {
344 t_coords(i) += t_n * t_ele_coords(i);
345 ++t_ele_coords;
346 ++t_n;
347 }
348
349 for (auto ii : {0, 1, 2})
350 t_coords(ii) = set_float_precision(t_coords(ii));
351
352 ++t_coords;
353 }
354
355 Tag th;
356 int def_in_the_loop = -1;
357 CHKERR postProcMesh.tag_get_handle("NB_IN_THE_LOOP", 1, MB_TYPE_INTEGER, th,
358 MB_TAG_CREAT | MB_TAG_SPARSE,
359 &def_in_the_loop);
360
361 postProcElements.clear();
362 const int num_el = level_ref_ele.size1();
363 const int num_nodes_on_ele = level_ref_ele.size2();
364 auto start_e = start_ele_handle + ele_count;
365 postProcElements = Range(start_e, start_e + num_el - 1);
366 for (auto tt = 0; tt != level_ref_ele.size1(); ++tt, ++ele_count) {
367 for (int nn = 0; nn != num_nodes_on_ele; ++nn) {
368 conn[num_nodes_on_ele * ele_count + nn] =
369 mapGaussPts[level_ref_ele(tt, nn)];
370 }
371 }
372
373 const int n_in_the_loop = E::nInTheLoop;
374 CHKERR postProcMesh.tag_clear_data(th, postProcElements, &n_in_the_loop);
376
378 };
379
380 CHKERR set_gauss_pts(
381
382 ref_ele->levelGaussPtsOnRefMesh, ref_ele->levelRef,
383 ref_ele->levelShapeFunctions,
384
385 ref_ele->startingVertEleHandle, ref_ele->startingEleHandle,
386 ref_ele->verticesOnEleArrays, ref_ele->eleConn, ref_ele->countVertEle,
387 ref_ele->countEle
388
389 );
390
392};
#define CHK_MOAB_THROW(err, msg)
Check error code of MoAB function and throw MoFEM exception.
Definition: definitions.h:576
FTensor::Index< 'i', SPACE_DIM > i
UBlasVector< double > VectorDouble
Definition: Types.hpp:68
auto type_from_handle(const EntityHandle h)
get type from entity handle
Definition: Templates.hpp:1594
virtual MoFEMErrorCode transferTags()

◆ setTagsToTransfer()

template<typename E >
MoFEMErrorCode MoFEM::PostProcBrokenMeshInMoabBase< E >::setTagsToTransfer ( std::vector< Tag >  tags_to_transfer)

Definition at line 172 of file PostProcBrokenMeshInMoabBase.hpp.

173 {
175 tagsToTransfer.swap(tags_to_transfer);
177}
std::vector< Tag > tagsToTransfer

◆ transferTags()

template<typename E >
MoFEMErrorCode MoFEM::PostProcBrokenMeshInMoabBase< E >::transferTags
protectedvirtual

transfer tags from mesh to post-process mesh

Definition at line 194 of file PostProcBrokenMeshInMoabBase.hpp.

194 {
196
197 auto &calc_mesh = this->mField.get_moab();
198
199 auto name = [&](auto tag) {
200 std::string name;
201 CHK_MOAB_THROW(calc_mesh.tag_get_name(tag, name), "get name");
202 return name;
203 };
204
205 auto data_type = [&](auto tag) {
206 moab::DataType data_type;
207 CHK_MOAB_THROW(calc_mesh.tag_get_data_type(tag, data_type), "get data type");
208 return data_type;
209 };
210
211 auto type = [&](auto tag) {
212 moab::TagType type;
213 CHK_MOAB_THROW(calc_mesh.tag_get_type(tag, type), "get tag type");
214 return type;
215 };
216
217 auto length = [&](auto tag) {
218 int length;
219 CHK_MOAB_THROW(calc_mesh.tag_get_length(tag, length), "get length ");
220 return length;
221 };
222
223 auto default_value = [&](auto tag) {
224 const void* def_val;
225 int size;
226 CHK_MOAB_THROW(calc_mesh.tag_get_default_value(tag, def_val, size),
227 "get default tag value");
228 return def_val;
229 };
230
231 auto data_ptr = [&](auto tag) {
232 const void *tag_data;
233 auto core_ent = this->getFEEntityHandle();
234 CHK_MOAB_THROW(calc_mesh.tag_get_by_ptr(tag, &core_ent, 1, &tag_data),
235 "get tag data");
236 return tag_data;
237 };
238
239 for (auto tag : tagsToTransfer) {
240 Tag tag_postproc;
241 CHKERR postProcMesh.tag_get_handle(
242 name(tag).c_str(), length(tag), data_type(tag), tag_postproc,
243 type(tag) | MB_TAG_CREAT, default_value(tag));
244 CHKERR postProcMesh.tag_clear_data(tag_postproc, postProcElements,
245 data_ptr(tag));
246 }
247
249};

Member Data Documentation

◆ coreMesh

template<typename E >
moab::Core MoFEM::PostProcBrokenMeshInMoabBase< E >::coreMesh
protected

Definition at line 159 of file PostProcBrokenMeshInMoabBase.hpp.

◆ mapGaussPts

template<typename E >
std::vector<EntityHandle> MoFEM::PostProcBrokenMeshInMoabBase< E >::mapGaussPts

Definition at line 100 of file PostProcBrokenMeshInMoabBase.hpp.

◆ postProcElements

template<typename E >
Range MoFEM::PostProcBrokenMeshInMoabBase< E >::postProcElements

Definition at line 101 of file PostProcBrokenMeshInMoabBase.hpp.

◆ postProcMesh

template<typename E >
moab::Interface& MoFEM::PostProcBrokenMeshInMoabBase< E >::postProcMesh
protected

Definition at line 160 of file PostProcBrokenMeshInMoabBase.hpp.

◆ refElementsMap

template<typename E >
std::map<EntityType, PostProcGenerateRefMeshPtr> MoFEM::PostProcBrokenMeshInMoabBase< E >::refElementsMap
protected

Definition at line 162 of file PostProcBrokenMeshInMoabBase.hpp.

◆ tagsToTransfer

template<typename E >
std::vector<Tag> MoFEM::PostProcBrokenMeshInMoabBase< E >::tagsToTransfer
protected

Set of tags on mesh to transfer to postprocessing mesh

Definition at line 164 of file PostProcBrokenMeshInMoabBase.hpp.


The documentation for this struct was generated from the following file: