v0.15.0
Loading...
Searching...
No Matches
Core.cpp
Go to the documentation of this file.
1/** \file Core.cpp
2 * \brief Multi-index containers, data structures and other low-level functions
3 */
4
5
6#include <MoFEM.hpp>
7
9
10extern "C" {
12}
13
14namespace MoFEM {
15
16WrapMPIComm::WrapMPIComm(MPI_Comm comm, bool petsc)
17 : comm(comm), isPetscComm(petsc) {
18 if (isPetscComm) {
19 ierr = PetscCommDuplicate(comm, &duplicatedComm, NULL);
20 CHKERRABORT(comm, ierr);
21 } else {
22 int ierr = MPI_Comm_dup(comm, &duplicatedComm);
23 if (ierr) {
24 THROW_MESSAGE("MPI_Comm_dup not working");
25 }
26 }
27}
29 if (isPetscComm) {
30 ierr = PetscCommDestroy(&duplicatedComm);
31 CHKERRABORT(comm, ierr);
32 } else {
33 int ierr = MPI_Comm_free(&duplicatedComm);
34 if (ierr) {
35 CHKERRABORT(comm, MOFEM_DATA_INCONSISTENCY);
36 }
37 }
38}
39
40constexpr const int CoreTmp<0>::value;
41constexpr const int CoreTmp<-1>::value;
42
43MoFEMErrorCode Core::query_interface(boost::typeindex::type_index type_index,
44 UnknownInterface **iface) const {
46 *iface = NULL;
47 if (type_index == boost::typeindex::type_id<CoreInterface>()) {
48 *iface = static_cast<CoreInterface *>(const_cast<Core *>(this));
50 } else if (type_index ==
51 boost::typeindex::type_id<DeprecatedCoreInterface>()) {
52 *iface = static_cast<DeprecatedCoreInterface *>(const_cast<Core *>(this));
54 }
55
56 // Get sub-interface
57 auto it = iFaces.find(type_index);
58 if (it != iFaces.end()) {
59 *iface = it->second;
61 }
62
63 *iface = NULL;
64 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown interface");
66}
67
70PetscBool Core::isInitialized;
71
72MoFEMErrorCode Core::Initialize(int *argc, char ***args, const char file[],
73 const char help[]) {
74
75 MPI_Initialized(&mpiInitialised);
76 if (!mpiInitialised)
77 MPI_Init(argc, args);
78
79 PetscInitialized(&isInitialized);
80 if (isInitialized == PETSC_FALSE) {
81 for (int i = 0; i < *argc - 1; ++i) {
82 if (std::strcmp((*args)[i], "-param_file") == 0) {
83 file = (*args)[i + 1];
84 break;
85 }
86 }
87 PetscInitialize(argc, args, file, help);
88 PetscPushErrorHandler(mofem_error_handler, PETSC_NULLPTR);
89 }
90
91 LogManager::createDefaultSinks(MPI_COMM_WORLD);
92 PetscVFPrintf = LogManager::logPetscFPrintf;
94 isGloballyInitialised = true;
95
96 MOFEM_LOG_CHANNEL("WORLD");
97 char petsc_version[255];
98 CHKERR PetscGetVersion(petsc_version, 255);
99 MOFEM_LOG_C("WORLD", Sev::inform, "MoFEM version %d.%d.%d (%s %s)",
100 MoFEM_VERSION_MAJOR, MoFEM_VERSION_MINOR, MoFEM_VERSION_BUILD,
101 MOAB_VERSION_STRING, petsc_version);
102 MOFEM_LOG_C("WORLD", Sev::inform, "git commit id %s", GIT_SHA1_NAME);
103
104 auto log_time = [&](const auto perefix, auto time) {
105 MOFEM_LOG("WORLD", Sev::inform)
106 << perefix << time.date().year() << "-" << time.date().month() << "-"
107 << time.date().day() << " " << time.time_of_day().hours() << ":"
108 << time.time_of_day().minutes() << ":" << time.time_of_day().seconds();
109 };
110
111 // Get current system time
112 log_time("Local time: ", boost::posix_time::second_clock::local_time());
113 log_time("UTC time: ", boost::posix_time::second_clock::universal_time());
114
115 return MOFEM_SUCCESS;
116}
117
119 if (isGloballyInitialised) {
120 PetscPopErrorHandler();
121 isGloballyInitialised = false;
122
123 if (isInitialized == PETSC_FALSE) {
124 PetscBool is_finalized;
125 PetscFinalized(&is_finalized);
126 if (!is_finalized)
127 PetscFinalize();
128 }
129
130 if (!mpiInitialised) {
131 int mpi_finalized;
132 MPI_Finalized(&mpi_finalized);
133 if (!mpi_finalized)
134 MPI_Finalize();
135 }
136 }
137
138 return 0;
139}
140
141// Use SFINAE to decide which template should be run,
142// if exist getSubInterfaceOptions run this one.
143template <class T>
144static auto get_sub_iface_options_imp(T *const ptr, int)
145 -> decltype(ptr->getSubInterfaceOptions()) {
146 return ptr->getSubInterfaceOptions();
147};
148
149// Use SFINAE to decide which template should be run,
150// if getSubInterfaceOptions not exist run this one.
151template <class T>
152static auto get_sub_iface_options_imp(T *const ptr, long) -> MoFEMErrorCode {
153 return 0;
154};
155
156template <class T>
157static auto get_event_options_imp(T *const ptr, int)
158 -> decltype(ptr->getEventOptions()) {
159 return ptr->getEventptions();
160};
161
162// Use SFINAE to decide which template should be run,
163// if getSubInterfaceOptions not exist run this one.
164// See SFINAE:
165// https://stackoverflow.com/questions/257288/is-it-possible-to-write-a-template-to-check-for-a-functions-existence
166// https://en.wikipedia.org/wiki/Substitution_failure_is_not_an_error
167template <class T>
168static auto get_event_options_imp(T *const ptr, long) -> MoFEMErrorCode {
169 return 0;
170};
171
172template <class IFACE> MoFEMErrorCode Core::regSubInterface() {
174 CHKERR registerInterface<IFACE>(true);
175 IFACE *ptr = new IFACE(*this);
176
177 // If sub interface has function getSubInterfaceOptions run
178 // it after construction. getSubInterfaceOptions is used to
179 // get parameters from command line.
180 auto get_sub_iface_options = [](auto *const ptr) {
181 return get_sub_iface_options_imp(ptr, 0);
182 };
183 CHKERR get_sub_iface_options(ptr);
184
185 auto type_idx = boost::typeindex::type_id<IFACE>();
186 iFaces.insert(type_idx, ptr);
188}
189
190template <class IFACE> MoFEMErrorCode Core::regEvents() {
192 auto ptr = boost::make_shared<IFACE>();
193 // See SFINAE:
194 // https://stackoverflow.com/questions/257288/is-it-possible-to-write-a-template-to-check-for-a-functions-existence
195 // https://en.wikipedia.org/wiki/Substitution_failure_is_not_an_error
196 auto get_event_options = [](auto *const ptr) {
197 return get_event_options_imp(ptr, 0);
198 };
199 CHKERR get_event_options(ptr.get());
201}
202
204 MPI_Comm comm, const int verbose) {
206
207 // This is deprecated ONE should use MoFEM::Core::Initialize
208 if (!isGloballyInitialised)
209 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
210 "MoFEM globally is not initialised, call MoFEM::Core::Initialize");
211
212 // Create duplicate communicator
213 wrapMPIMOABComm = boost::make_shared<WrapMPIComm>(comm, false);
214
215 MPI_Comm_size(mofemComm, &sIze);
216 MPI_Comm_rank(mofemComm, &rAnk);
217
218 // CHeck if moab has set communicator if not set communicator internally
219 ParallelComm *pComm = ParallelComm::get_pcomm(&moab, MYPCOMM_INDEX);
220 if (pComm == NULL)
221 pComm = new ParallelComm(&moab, wrapMPIMOABComm->get_comm());
222
223 // Register interfaces for this implementation
224 CHKERR registerInterface<UnknownInterface>();
225 CHKERR registerInterface<CoreInterface>();
226 CHKERR registerInterface<DeprecatedCoreInterface>();
227
228 // Register MOFEM events in PETSc
229 PetscLogEventRegister("FE_preProcess", 0, &MOFEM_EVENT_preProcess);
230 PetscLogEventRegister("FE_operator", 0, &MOFEM_EVENT_operator);
231 PetscLogEventRegister("FE_postProcess", 0, &MOFEM_EVENT_postProcess);
232 PetscLogEventRegister("MoFEMCreateMat", 0, &MOFEM_EVENT_createMat);
233
234 MOFEM_LOG_CHANNEL("WORLD");
235 MOFEM_LOG_CHANNEL("SELF");
236 MOFEM_LOG_CHANNEL("SYNC");
237
239}
240
241Core::CoreTmp(moab::Interface &moab, ///< MoAB interface
242 MPI_Comm comm, ///< MPI communicator
243 const int verbose ///< Verbosity level
244
245 )
246 : CoreTmp(moab, comm, verbose, CoreValue<0>()) {
247
248 // Register sub-interfaces
249 ierr = this->registerSubInterfaces();
250 CHKERRABORT(comm, ierr);
251 ierr = this->clearMap();
252 CHKERRABORT(comm, ierr);
253 ierr = this->getTags();
254 CHKERRABORT(comm, ierr);
255 ierr = this->getOptions(verbose);
256 CHKERRABORT(comm, ierr);
257
258 this->basicEntityDataPtr = boost::make_shared<BasicEntityData>(moab);
259 setRefEntBasicDataPtr(*this, this->basicEntityDataPtr);
260
261 ierr = this->initialiseDatabaseFromMesh(verbose);
262 CHKERRABORT(comm, ierr);
263}
264
265CoreTmp<-1>::CoreTmp(moab::Interface &moab, ///< MoAB interface
266 MPI_Comm comm, ///< MPI communicator
267 const int verbose ///< Verbosity level
268
269 )
270 : CoreTmp<0>(moab, comm, verbose, CoreValue<-1>()) {
271
272 // Register sub-interfaces
273 ierr = this->registerSubInterfaces();
274 CHKERRABORT(comm, ierr);
275 ierr = this->clearMap();
276 CHKERRABORT(comm, ierr);
277 ierr = this->getTags();
278 CHKERRABORT(comm, ierr);
279 ierr = this->getOptions(verbose);
280 CHKERRABORT(comm, ierr);
281
282 this->basicEntityDataPtr = boost::make_shared<BasicEntityData>(moab);
283 setRefEntBasicDataPtr(*this, this->basicEntityDataPtr);
284
285 ierr = this->initialiseDatabaseFromMesh(verbose);
286 CHKERRABORT(comm, ierr);
287}
288
290 PetscBool is_finalized;
291 PetscFinalized(&is_finalized);
292 // Destroy interfaces
293 iFaces.clear();
294 // This is deprecated ONE should use MoFEM::Core::Initialize
295 if (isGloballyInitialised && is_finalized) {
296 isGloballyInitialised = false;
297 }
298}
299
301 BOOST_LOG_SCOPED_THREAD_ATTR("Timeline", attrs::timer());
302 MOFEM_LOG_CHANNEL("WORLD");
303
305 if (verb == -1)
306 verb = verbose;
307
308 Range ref_elems_to_add;
309
310 MOFEM_LOG("WORLD", Sev::verbose) << "Get MoFEM meshsets";
311 // Initialize database
312 Range meshsets;
313 CHKERR get_moab().get_entities_by_type(0, MBENTITYSET, meshsets, false);
314 Range special_meshsets;
315 for (auto mit : meshsets) {
316 BitFieldId field_id;
317 // Get bit id form field tag
318 CHKERR get_moab().tag_get_data(th_FieldId, &mit, 1, &field_id);
319 // Check if meshset if field meshset
320 if (field_id != 0) {
321
322 const void *tag_name;
323 int tag_name_size;
324 CHKERR get_moab().tag_get_by_ptr(
325 th_FieldName, &mit, 1, (const void **)&tag_name, &tag_name_size);
326
327 if (verb > QUIET)
328 MOFEM_LOG("WORLD", Sev::verbose)
329 << "Read field "
330 << boost::string_ref((char *)tag_name, tag_name_size);
331
332 auto p = fIelds.insert(boost::make_shared<Field>(moab, mit));
333
334 if (!p.second) {
335 // Field meshset exists, remove duplicate meshsets from other
336 // processors.
337 Range ents;
338 CHKERR get_moab().get_entities_by_handle(mit, ents, true);
339 CHKERR get_moab().add_entities((*p.first)->getMeshset(), ents);
340 CHKERR get_moab().delete_entities(&mit, 1);
341 } else {
342 special_meshsets.insert(mit);
343 }
344 }
345 // Check for finite elements
346 BitFieldId fe_id;
347 // Get bit id from fe tag
348 CHKERR get_moab().tag_get_data(th_FEId, &mit, 1, &fe_id);
349 // check if meshset is finite element meshset
350 if (fe_id != 0) {
351 std::pair<FiniteElement_multiIndex::iterator, bool> p =
352 finiteElements.insert(
353 boost::shared_ptr<FiniteElement>(new FiniteElement(moab, mit)));
354 if (verb > QUIET)
355 MOFEM_LOG("WORLD", Sev::verbose) << "Read finite element " << **p.first;
356
357 Range ents;
358 CHKERR get_moab().get_entities_by_type(mit, MBENTITYSET, ents, false);
359 CHKERR get_moab().get_entities_by_handle(mit, ents, true);
360 ref_elems_to_add.merge(ents);
361 if (!p.second) {
362 // Finite element mesh set exist, could be created on other processor.
363 // Remove duplicate.
364 CHKERR get_moab().add_entities((*p.first)->getMeshset(), ents);
365 CHKERR get_moab().delete_entities(&mit, 1);
366 } else {
367 special_meshsets.insert(mit);
368 }
369 }
370 BitProblemId problem_id;
371 // get bit id form problem tag
372 CHKERR get_moab().tag_get_data(th_ProblemId, &mit, 1, &problem_id);
373 // check if meshset if problem meshset
374 if (problem_id != 0) {
375 std::pair<Problem_multiIndex::iterator, bool> p =
376 pRoblems.insert(Problem(moab, mit));
377 if (verb > QUIET) {
378 MOFEM_LOG("WORLD", Sev::verbose) << "Read problem " << *p.first;
379 MOFEM_LOG("WORLD", Sev::noisy)
380 << "\tBitRef " << p.first->getBitRefLevel() << " BitMask "
381 << p.first->getBitRefLevelMask();
382 }
383
384 if (!p.second) {
385 // Problem meshset exists, could be created on other processor.
386 // Remove duplicate.
387 Range ents;
388 CHKERR get_moab().get_entities_by_handle(mit, ents, true);
389 CHKERR get_moab().get_entities_by_type(mit, MBENTITYSET, ents, true);
390 CHKERR get_moab().add_entities(p.first->meshset, ents);
391 CHKERR get_moab().delete_entities(&mit, 1);
392 } else {
393 special_meshsets.insert(mit);
394 }
395 }
396 }
397 MOFEM_LOG("WORLD", Sev::verbose) << "Get MoFEM meshsets <- done";
398
399 // Add entities to database
400 MOFEM_LOG("WORLD", Sev::verbose) << "Add entities to database";
401 Range bit_ref_ents;
402 CHKERR get_moab().get_entities_by_handle(0, bit_ref_ents, false);
403 bit_ref_ents = subtract(bit_ref_ents, special_meshsets);
404 CHKERR getInterface<BitRefManager>()->filterEntitiesByRefLevel(
405 BitRefLevel().set(), BitRefLevel().set(), bit_ref_ents);
406 CHKERR getInterface<BitRefManager>()->setEntitiesBitRefLevel(bit_ref_ents);
407 CHKERR getInterface<BitRefManager>()->setElementsBitRefLevel(
408 ref_elems_to_add);
409 MOFEM_LOG("WORLD", Sev::verbose) << "Add entities to database <- done";
410
411 // Build field entities
412 MOFEM_LOG("WORLD", Sev::verbose) << "Add field to database";
413 for (auto field : fIelds) {
414 if (field->getSpace() != NOSPACE) {
415 Range ents_of_id_meshset;
416 CHKERR get_moab().get_entities_by_handle(field->getMeshset(),
417 ents_of_id_meshset, false);
418 CHKERR set_field_order(ents_of_id_meshset, field->getId(), -1, verb);
419 }
420 }
421 MOFEM_LOG("WORLD", Sev::verbose) << "Add field to database <- done";
422
423 if (initaliseAndBuildField || initaliseAndBuildFiniteElements) {
424 MOFEM_LOG("WORLD", Sev::verbose) << "Build fields elements";
425 CHKERR build_fields(verb);
426 MOFEM_LOG("WORLD", Sev::verbose) << "Build fields elements <- done";
427 if (initaliseAndBuildFiniteElements) {
428 MOFEM_LOG("WORLD", Sev::verbose) << "Build finite elements";
429 CHKERR build_finite_elements(verb);
430 MOFEM_LOG("WORLD", Sev::verbose) << "Build finite elements <- done";
431 }
432 }
433
434 if (verb > VERY_NOISY) {
435 list_fields();
436 list_finite_elements();
437 list_problem();
438 }
439
440 // Initialize interfaces
441 MOFEM_LOG("WORLD", Sev::verbose) << "Initialise interfaces from mesh";
442
443 MOFEM_LOG("WORLD", Sev::verbose) << "Initialise MeshsetManager";
444 CHKERR getInterface<MeshsetsManager>() -> initialiseDatabaseFromMesh(verb);
445 MOFEM_LOG("WORLD", Sev::verbose) << "Initialise MeshsetManager <- done";
446 MOFEM_LOG("WORLD", Sev::verbose) << "Initialise SeriesRecorder";
447 CHKERR getInterface<SeriesRecorder>() -> initialiseDatabaseFromMesh(verb);
448 MOFEM_LOG("WORLD", Sev::verbose) << "Initialise SeriesRecorder <- done";
449
450 MOFEM_LOG("WORLD", Sev::verbose) << "Initialise interfaces from mesh <- done";
451
453}
454
455MoFEMErrorCode Core::setMoabInterface(moab::Interface &new_moab, int verb) {
457 if (verb == -1)
458 verb = verbose;
459
460 // clear moab database
461 CHKERR clearMap();
462
463 // set new reference
464 moab = std::ref(new_moab);
465
466 // check if moab has set communicator if not set communicator internally
467 ParallelComm *pComm = ParallelComm::get_pcomm(&new_moab, MYPCOMM_INDEX);
468 if (pComm == NULL) {
469 pComm = new ParallelComm(&new_moab, wrapMPIMOABComm->get_comm());
470 }
471
472 // create MoFEM tags
473 CHKERR getTags();
474
475 // Create basic entity data struture
476 basicEntityDataPtr = boost::make_shared<BasicEntityData>(moab);
477 setRefEntBasicDataPtr(*this, this->basicEntityDataPtr);
478
479 // Initalise database
480 CHKERR this->initialiseDatabaseFromMesh(verb);
481
483};
484
487
488 iFaces.clear();
489
490 // Register sub interfaces
491 CHKERR regSubInterface<LogManager>();
492 CHKERR regSubInterface<Simple>();
493 CHKERR regSubInterface<OperatorsTester>();
494 CHKERR regSubInterface<PipelineManager>();
495 CHKERR regSubInterface<ProblemsManager>();
496 CHKERR regSubInterface<MatrixManager>();
497 CHKERR regSubInterface<ISManager>();
498 CHKERR regSubInterface<VecManager>();
499 CHKERR regSubInterface<FieldBlas>();
500 CHKERR regSubInterface<BitRefManager>();
501 CHKERR regSubInterface<Tools>();
502 CHKERR regSubInterface<CommInterface>();
503 CHKERR regSubInterface<MeshsetsManager>();
504 CHKERR regSubInterface<NodeMergerInterface>();
505 CHKERR regSubInterface<PrismsFromSurfaceInterface>();
506 CHKERR regSubInterface<MeshRefinement>();
507 CHKERR regSubInterface<PrismInterface>();
508 CHKERR regSubInterface<CutMeshInterface>();
509 CHKERR regSubInterface<SeriesRecorder>();
510#ifdef WITH_TETGEN
511 CHKERR regSubInterface<TetGenInterface>();
512#endif
513#ifdef WITH_MED
514 CHKERR regSubInterface<MedInterface>();
515#endif
516 CHKERR regSubInterface<FieldEvaluatorInterface>();
517 CHKERR regSubInterface<BcManager>();
518
519 // Register events
520 CHKERR regEvents<SchurEvents>();
521
523};
524
527 // Cleaning databases in interfaces
528 CHKERR getInterface<SeriesRecorder>()->clearMap();
529 CHKERR getInterface<MeshsetsManager>()->clearMap();
530 CHKERR getInterface<CutMeshInterface>()->clearMap();
531 // Cleaning databases
532 refinedEntities.clear();
533 refinedFiniteElements.clear();
534 fIelds.clear();
535 entsFields.clear();
536 dofsField.clear();
537 finiteElements.clear();
538 entsFiniteElements.clear();
539 entFEAdjacencies.clear();
540 pRoblems.clear();
542}
543
546 if (verb == -1)
547 verb = verbose;
548 std::pair<RefEntity_multiIndex::iterator, bool> p_ent;
549 p_ent = refinedEntities.insert(
550 boost::make_shared<RefEntity>(basicEntityDataPtr, prism));
551 if (p_ent.second) {
552 std::pair<RefElement_multiIndex::iterator, bool> p;
553 p = refinedFiniteElements.insert(
554 boost::shared_ptr<RefElement>(new RefElement_PRISM(*p_ent.first)));
555 int num_nodes;
556 const EntityHandle *conn;
557 CHKERR get_moab().get_connectivity(prism, conn, num_nodes, true);
558 Range face_side3, face_side4;
559 CHKERR get_moab().get_adjacencies(conn, 3, 2, false, face_side3);
560 CHKERR get_moab().get_adjacencies(&conn[3], 3, 2, false, face_side4);
561 if (face_side3.size() != 1)
562 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
563 "prism don't have side face 3");
564 if (face_side4.size() != 1)
565 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
566 "prims don't have side face 4");
567 p.first->get()->getSideNumberPtr(*face_side3.begin());
568 p.first->get()->getSideNumberPtr(*face_side4.begin());
569 }
571}
572
575
576 const EntityHandle root_meshset = get_moab().get_root_set();
577 if (root_meshset) {
578 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
579 "Root meshset should be 0");
580 }
581
582 // Set version
583 {
584 Version version;
585 CHKERR getFileVersion(moab, version);
586 }
587
588 // Global Variables
589 {
590
591 auto check_tag_allocated = [](auto &rval) {
593 if (rval == MB_ALREADY_ALLOCATED)
594 rval = MB_SUCCESS;
595 else
596 CHKERRG(rval);
598 };
599
600 // Safety nets
601 int def_bool = 0;
602 rval = get_moab().tag_get_handle("_MoFEMBuild", 1, MB_TYPE_INTEGER,
603 th_MoFEMBuild, MB_TAG_CREAT | MB_TAG_MESH,
604 &def_bool);
605 CHKERR check_tag_allocated(rval);
606
607 CHKERR get_moab().tag_get_by_ptr(th_MoFEMBuild, &root_meshset, 1,
608 (const void **)&buildMoFEM);
609 }
610
611 // Tags saved in vtk-files
612 {
613 const int def_part = -1;
614 CHKERR get_moab().tag_get_handle("PARTITION", 1, MB_TYPE_INTEGER, th_Part,
615 MB_TAG_CREAT | MB_TAG_SPARSE, &def_part);
616 }
617
618 // Tags Ref
619 {
620
621 // Fix size of bir ref level tags
623
624 const int def_part = -1;
625 CHKERR get_moab().tag_get_handle("_MeshsetPartition", 1, MB_TYPE_INTEGER,
626 th_Part, MB_TAG_CREAT | MB_TAG_SPARSE,
627 &def_part);
628 EntityHandle def_handle = 0;
629 CHKERR get_moab().tag_get_handle("_RefParentHandle", 1, MB_TYPE_HANDLE,
630 th_RefParentHandle,
631 MB_TAG_CREAT | MB_TAG_SPARSE, &def_handle);
632 BitRefLevel def_bit_level = 0;
633 CHKERR get_moab().tag_get_handle(
634 "_RefBitLevel", sizeof(BitRefLevel), MB_TYPE_OPAQUE, th_RefBitLevel,
635 MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_SPARSE, &def_bit_level);
636 BitRefLevel def_bit_level_mask = BitRefLevel().set();
637 CHKERR get_moab().tag_get_handle(
638 "_RefBitLevelMask", sizeof(BitRefLevel), MB_TYPE_OPAQUE,
639 th_RefBitLevel_Mask, MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_SPARSE,
640 &def_bit_level_mask);
641 BitRefEdges def_bit_edge = 0;
642 CHKERR get_moab().tag_get_handle(
643 "_RefBitEdge", sizeof(BitRefEdges), MB_TYPE_OPAQUE, th_RefBitEdge,
644 MB_TAG_CREAT | MB_TAG_SPARSE | MB_TAG_BYTES, &def_bit_edge);
645 }
646
647 // Tags Field
648 {
649 const unsigned long int def_id = 0;
650 CHKERR get_moab().tag_get_handle(
651 "_FieldId", sizeof(BitFieldId), MB_TYPE_OPAQUE, th_FieldId,
652 MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_SPARSE, &def_id);
653 FieldSpace def_space = LASTSPACE;
654 CHKERR get_moab().tag_get_handle(
655 "_FieldSpace", sizeof(FieldSpace), MB_TYPE_OPAQUE, th_FieldSpace,
656 MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_SPARSE, &def_space);
657 FieldContinuity def_continuity = LASTCONTINUITY;
658 CHKERR get_moab().tag_get_handle(
659 "_FieldContinuity", sizeof(FieldContinuity), MB_TYPE_OPAQUE,
660 th_FieldContinuity, MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_SPARSE,
661 &def_continuity);
663 CHKERR get_moab().tag_get_handle(
664 "_FieldBase", sizeof(FieldApproximationBase), MB_TYPE_OPAQUE,
665 th_FieldBase, MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_SPARSE, &def_base);
666 const int def_val_len = 0;
667 CHKERR get_moab().tag_get_handle(
668 "_FieldName", def_val_len, MB_TYPE_OPAQUE, th_FieldName,
669 MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_VARLEN | MB_TAG_SPARSE, NULL);
670 CHKERR get_moab().tag_get_handle(
671 "_FieldName_DataNamePrefix", def_val_len, MB_TYPE_OPAQUE,
672 th_FieldName_DataNamePrefix,
673 MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_VARLEN | MB_TAG_SPARSE, NULL);
674 }
675
676 // Tags FE
677 {
678 const unsigned long int def_id = 0;
679 const int def_val_len = 0;
680 CHKERR get_moab().tag_get_handle(
681 "_FEId", sizeof(BitFEId), MB_TYPE_OPAQUE, th_FEId,
682 MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_SPARSE, &def_id);
683 CHKERR get_moab().tag_get_handle(
684 "_FEName", def_val_len, MB_TYPE_OPAQUE, th_FEName,
685 MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_VARLEN | MB_TAG_SPARSE, NULL);
686 CHKERR get_moab().tag_get_handle(
687 "_FEIdCol", sizeof(BitFieldId), MB_TYPE_OPAQUE, th_FEIdCol,
688 MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_SPARSE, &def_id);
689 CHKERR get_moab().tag_get_handle(
690 "_FEIdRow", sizeof(BitFieldId), MB_TYPE_OPAQUE, th_FEIdRow,
691 MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_SPARSE, &def_id);
692 CHKERR get_moab().tag_get_handle(
693 "_FEIdData", sizeof(BitFieldId), MB_TYPE_OPAQUE, th_FEIdData,
694 MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_SPARSE, &def_id);
695 }
696
697 // Tags Problem
698 {
699 const unsigned long int def_id = 0;
700 const int def_val_len = 0;
701 CHKERR get_moab().tag_get_handle(
702 "_ProblemId", sizeof(BitProblemId), MB_TYPE_OPAQUE, th_ProblemId,
703 MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_SPARSE, &def_id);
704 CHKERR get_moab().tag_get_handle(
705 "_ProblemFEId", sizeof(BitFEId), MB_TYPE_OPAQUE, th_ProblemFEId,
706 MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_SPARSE, &def_id);
707 CHKERR get_moab().tag_get_handle(
708 "_ProblemName", def_val_len, MB_TYPE_OPAQUE, th_ProblemName,
709 MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_VARLEN | MB_TAG_SPARSE, NULL);
710 DofIdx def_nbdofs = 0;
711 CHKERR get_moab().tag_get_handle(
712 "_ProblemNbDofsRow", sizeof(DofIdx), MB_TYPE_OPAQUE,
713 th_ProblemNbDofsRow, MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_SPARSE,
714 &def_nbdofs);
715 CHKERR get_moab().tag_get_handle(
716 "_ProblemNbDofsCol", sizeof(DofIdx), MB_TYPE_OPAQUE,
717 th_ProblemNbDofsCol, MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_SPARSE,
718 &def_nbdofs);
719 CHKERR get_moab().tag_get_handle(
720 "_ProblemLocalNbDofsRow", sizeof(DofIdx), MB_TYPE_OPAQUE,
721 th_ProblemLocalNbDofRow, MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_SPARSE,
722 &def_nbdofs);
723 CHKERR get_moab().tag_get_handle(
724 "_ProblemGhostNbDofsRow", sizeof(DofIdx), MB_TYPE_OPAQUE,
725 th_ProblemGhostNbDofRow, MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_SPARSE,
726 &def_nbdofs);
727 CHKERR get_moab().tag_get_handle(
728 "_ProblemLocalNbDofsCol", sizeof(DofIdx), MB_TYPE_OPAQUE,
729 th_ProblemLocalNbDofCol, MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_SPARSE,
730 &def_nbdofs);
731 CHKERR get_moab().tag_get_handle(
732 "_ProblemGhostNbDofsCol", sizeof(DofIdx), MB_TYPE_OPAQUE,
733 th_ProblemGhostNbDofCol, MB_TAG_CREAT | MB_TAG_BYTES | MB_TAG_SPARSE,
734 &def_nbdofs);
735 }
736
737 // Meshsets with boundary conditions and material sets
738 MeshsetsManager *meshsets_manager_ptr;
739 CHKERR getInterface(meshsets_manager_ptr);
740 CHKERR meshsets_manager_ptr->getTags(verb);
741
742 // Series recorder
743 SeriesRecorder *series_recorder_ptr;
744 CHKERR getInterface(series_recorder_ptr);
745 CHKERR series_recorder_ptr->getTags(verb);
746
748}
749
752 if (verb == -1)
753 verb = verbose;
754 CHKERR clearMap();
756}
757
760 if (verb == -1)
761 verb = verbose;
762 CHKERR this->clearMap();
763 CHKERR this->getTags(verb);
764 CHKERR this->initialiseDatabaseFromMesh(verb);
766}
767
768MoFEMErrorCode Core::set_moab_interface(moab::Interface &new_moab, int verb) {
769 return this->setMoabInterface(new_moab, verb);
770};
771
773 int verb) {
774 return this->setMoabInterface(new_moab, verb);
775};
776
779 if (verb == -1)
780 verb = verbose;
781
782 PetscOptionsBegin(mofemComm, optionsPrefix.c_str(), "Mesh cut options",
783 "See MoFEM documentation");
784
785 CHKERR PetscOptionsBool(
786 "-mofem_init_fields", "Initialise fields on construction", "",
787 initaliseAndBuildField, &initaliseAndBuildField, NULL);
788
789 CHKERR PetscOptionsBool(
790 "-mofem_init_fields", "Initialise fields on construction", "",
791 initaliseAndBuildFiniteElements, &initaliseAndBuildFiniteElements, NULL);
792
793 // TODO: Add read verbosity level
794 // TODO: Add option to initalise problems ??? - DO WE REALLY NEED THAT
795
796 PetscOptionsEnd();
797
799}
800
801// cubit meshsets
802
805 *fields_ptr = &fIelds;
807}
808
810Core::get_ref_ents(const RefEntity_multiIndex **refined_entities_ptr) const {
812 *refined_entities_ptr = &refinedEntities;
814}
816 const RefElement_multiIndex **refined_finite_elements_ptr) const {
818 *refined_finite_elements_ptr = &refinedFiniteElements;
820}
821
822MoFEMErrorCode Core::get_problem(const std::string &problem_name,
823 const Problem **problem_ptr) const {
825 typedef Problem_multiIndex::index<Problem_mi_tag>::type ProblemsByName;
826 const ProblemsByName &problems = pRoblems.get<Problem_mi_tag>();
827 ProblemsByName::iterator p_miit = problems.find(problem_name);
828 if (p_miit == problems.end()) {
829 SETERRQ(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
830 "problem < %s > not found, (top tip: check spelling)",
831 problem_name.c_str());
832 }
833 *problem_ptr = &*p_miit;
835}
836
838Core::get_problems(const Problem_multiIndex **problems_ptr) const {
840 *problems_ptr = &pRoblems;
842}
843
847 *field_ents = &entsFields;
849}
852 *dofs_ptr = &dofsField;
854}
855
859 *fe_ptr = &finiteElements;
861}
862
864 const EntFiniteElement_multiIndex **fe_ent_ptr) const {
866 *fe_ent_ptr = &entsFiniteElements;
868}
869
871 MeshsetsManager *meshsets_manager_ptr;
872 getInterface(meshsets_manager_ptr);
873 return meshsets_manager_ptr;
874}
875
877 MeshsetsManager *meshsets_manager_ptr;
878 getInterface(meshsets_manager_ptr);
879 return meshsets_manager_ptr;
880}
881
884 *dofs_elements_adjacency) const {
886 *dofs_elements_adjacency = &entFEAdjacencies;
888}
889
892 return &entFEAdjacencies;
893}
894
895const Field_multiIndex *Core::get_fields() const { return &fIelds; }
897 return &refinedEntities;
898}
900 return &refinedFiniteElements;
901}
903 return &finiteElements;
904}
906 return &entsFiniteElements;
907}
909 return &entsFields;
910}
911const DofEntity_multiIndex *Core::get_dofs() const { return &dofsField; }
912const Problem *Core::get_problem(const std::string problem_name) const {
913 const Problem *prb;
914 CHK_THROW_MESSAGE(get_problem(problem_name, &prb),
915 "Problem of given name not found");
916 return prb;
917}
918const Problem_multiIndex *Core::get_problems() const { return &pRoblems; }
919
920template <int V, typename std::enable_if<(V >= 0), int>::type * = nullptr>
921void set_ref_ent_basic_data_ptr_impl(boost::shared_ptr<BasicEntityData> &ptr) {
923};
924
925template <int V, typename std::enable_if<(V < 0), int>::type * = nullptr>
926void set_ref_ent_basic_data_ptr_impl(boost::shared_ptr<BasicEntityData> &ptr) {
927 return;
928};
929
930void Core::setRefEntBasicDataPtr(MoFEM::Interface &m_field,
931 boost::shared_ptr<BasicEntityData> &ptr) {
932
933 switch (m_field.getValue()) {
934 case -1:
935 set_ref_ent_basic_data_ptr_impl<-1>(ptr);
936 break;
937 case 0:
938 set_ref_ent_basic_data_ptr_impl<0>(ptr);
939 break;
940 case 1:
941 set_ref_ent_basic_data_ptr_impl<1>(ptr);
942 break;
943 default:
944 THROW_MESSAGE("Core index can vary from -1 to MAX_CORE_TMP");
945 }
946};
947
948boost::shared_ptr<RefEntityTmp<0>>
949Core::makeSharedRefEntity(MoFEM::Interface &m_field, const EntityHandle ent) {
950
951 boost::shared_ptr<RefEntityTmp<0>> ref_ent_ptr;
952
953 switch (m_field.getValue()) {
954 case -1:
955 ref_ent_ptr = boost::shared_ptr<RefEntityTmp<0>>(
956
957 new RefEntityTmp<-1>(m_field.get_basic_entity_data_ptr(), ent)
958
959 );
960 break;
961 case 0:
962 ref_ent_ptr = boost::shared_ptr<RefEntityTmp<0>>(
963
964 new RefEntityTmp<0>(m_field.get_basic_entity_data_ptr(), ent)
965
966 );
967 break;
968 case 1:
969 ref_ent_ptr = boost::shared_ptr<RefEntityTmp<0>>(
970
971 new RefEntityTmp<1>(m_field.get_basic_entity_data_ptr(), ent)
972
973 );
974 break;
975 default:
976 THROW_MESSAGE("Core index can vary from -1 to MAX_CORE_TMP");
977 }
978
979 return ref_ent_ptr;
980}
981
982boost::shared_ptr<RefEntityTmp<0>>
983Core::make_shared_ref_entity(const EntityHandle ent) {
984 return this->makeSharedRefEntity(*this, ent);
985}
986
987boost::shared_ptr<RefEntityTmp<0>>
989 return this->makeSharedRefEntity(*this, ent);
990}
991
992} // namespace MoFEM
multi_index_container< FieldEntityEntFiniteElementAdjacencyMap, indexed_by< ordered_unique< tag< Composite_Unique_mi_tag >, composite_key< FieldEntityEntFiniteElementAdjacencyMap, const_mem_fun< FieldEntityEntFiniteElementAdjacencyMap, UId, &FieldEntityEntFiniteElementAdjacencyMap::getEntUniqueId >, const_mem_fun< FieldEntityEntFiniteElementAdjacencyMap, UId, &FieldEntityEntFiniteElementAdjacencyMap::getFeUniqueId > > >, ordered_non_unique< tag< Unique_mi_tag >, const_mem_fun< FieldEntityEntFiniteElementAdjacencyMap, UId, &FieldEntityEntFiniteElementAdjacencyMap::getEntUniqueId > >, ordered_non_unique< tag< FE_Unique_mi_tag >, const_mem_fun< FieldEntityEntFiniteElementAdjacencyMap, UId, &FieldEntityEntFiniteElementAdjacencyMap::getFeUniqueId > >, ordered_non_unique< tag< FEEnt_mi_tag >, const_mem_fun< FieldEntityEntFiniteElementAdjacencyMap, EntityHandle, &FieldEntityEntFiniteElementAdjacencyMap::getFeHandle > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< FieldEntityEntFiniteElementAdjacencyMap, EntityHandle, &FieldEntityEntFiniteElementAdjacencyMap::getEntHandle > > > > FieldEntityEntFiniteElementAdjacencyMap_multiIndex
MultiIndex container keeps Adjacencies Element and dof entities adjacencies and vice versa.
void macro_is_deprecated_using_deprecated_function()
Is used to indicate that macro is deprecated Do nothing just triggers error at the compilation.
Definition Core.cpp:11
static PetscErrorCode mofem_error_handler(MPI_Comm comm, int line, const char *fun, const char *file, PetscErrorCode n, PetscErrorType p, const char *mess, void *ctx)
multi_index_container< boost::shared_ptr< Field >, indexed_by< hashed_unique< tag< BitFieldId_mi_tag >, const_mem_fun< Field, const BitFieldId &, &Field::getId >, HashBit< BitFieldId >, EqBit< BitFieldId > >, ordered_unique< tag< Meshset_mi_tag >, member< Field, EntityHandle, &Field::meshSet > >, ordered_unique< tag< FieldName_mi_tag >, const_mem_fun< Field, boost::string_ref, &Field::getNameRef > >, ordered_non_unique< tag< BitFieldId_space_mi_tag >, const_mem_fun< Field, FieldSpace, &Field::getSpace > > > > Field_multiIndex
Field_multiIndex for Field.
#define MOFEM_LOG_C(channel, severity, format,...)
static char help[]
@ QUIET
@ VERY_NOISY
FieldApproximationBase
approximation base
Definition definitions.h:58
@ LASTBASE
Definition definitions.h:69
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
FieldSpace
approximation spaces
Definition definitions.h:82
@ LASTSPACE
FieldSpace in [ 0, LASTSPACE )
Definition definitions.h:89
@ NOSPACE
Definition definitions.h:83
#define MYPCOMM_INDEX
default communicator number PCOMM
FieldContinuity
Field continuity.
Definition definitions.h:99
@ LASTCONTINUITY
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
@ MOFEM_OPERATION_UNSUCCESSFUL
Definition definitions.h:34
@ MOFEM_DATA_INCONSISTENCY
Definition definitions.h:31
@ MOFEM_SUCCESS
Definition definitions.h:30
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
multi_index_container< boost::shared_ptr< DofEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, const_mem_fun< DofEntity, UId, &DofEntity::getLocalUniqueId > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< DofEntity, EntityHandle, &DofEntity::getEnt > > > > DofEntity_multiIndex
MultiIndex container keeps DofEntity.
multi_index_container< boost::shared_ptr< EntFiniteElement >, indexed_by< ordered_unique< tag< Unique_mi_tag >, const_mem_fun< EntFiniteElement, UId, &EntFiniteElement::getLocalUniqueId > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< EntFiniteElement::interface_type_RefEntity, EntityHandle, &EntFiniteElement::getEnt > > > > EntFiniteElement_multiIndex
MultiIndex container for EntFiniteElement.
multi_index_container< Problem, indexed_by< ordered_unique< tag< Meshset_mi_tag >, member< Problem, EntityHandle, &Problem::meshset > >, hashed_unique< tag< BitProblemId_mi_tag >, const_mem_fun< Problem, BitProblemId, &Problem::getId >, HashBit< BitProblemId >, EqBit< BitProblemId > >, hashed_unique< tag< Problem_mi_tag >, const_mem_fun< Problem, std::string, &Problem::getName > > > > Problem_multiIndex
MultiIndex for entities for Problem.
multi_index_container< boost::shared_ptr< FiniteElement >, indexed_by< hashed_unique< tag< FiniteElement_Meshset_mi_tag >, member< FiniteElement, EntityHandle, &FiniteElement::meshset > >, hashed_unique< tag< BitFEId_mi_tag >, const_mem_fun< FiniteElement, BitFEId, &FiniteElement::getId >, HashBit< BitFEId >, EqBit< BitFEId > >, ordered_unique< tag< FiniteElement_name_mi_tag >, const_mem_fun< FiniteElement, boost::string_ref, &FiniteElement::getNameRef > > > > FiniteElement_multiIndex
MultiIndex for entities for FiniteElement.
#define MOFEM_LOG(channel, severity)
Log.
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
MoFEMErrorCode getTags(int verb=-1)
Get tag handles used on meshsets.
FTensor::Index< 'i', SPACE_DIM > i
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
std::bitset< BITFEID_SIZE > BitFEId
Finite element Id.
Definition Types.hpp:43
int DofIdx
Index of DOF.
Definition Types.hpp:18
std::bitset< BITPROBLEMID_SIZE > BitProblemId
Problem Id.
Definition Types.hpp:44
std::bitset< BITFIELDID_SIZE > BitFieldId
Field Id.
Definition Types.hpp:42
std::bitset< BITREFEDGES_SIZE > BitRefEdges
Definition Types.hpp:34
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition Types.hpp:40
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
multi_index_container< boost::shared_ptr< RefEntity >, indexed_by< ordered_unique< tag< Ent_mi_tag >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getEnt > >, ordered_non_unique< tag< Ent_Ent_mi_tag >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > >, ordered_non_unique< tag< Composite_EntType_and_ParentEntType_mi_tag >, composite_key< RefEntity, const_mem_fun< RefEntity, EntityType, &RefEntity::getEntType >, const_mem_fun< RefEntity, EntityType, &RefEntity::getParentEntType > > >, ordered_non_unique< tag< Composite_ParentEnt_And_EntType_mi_tag >, composite_key< RefEntity, const_mem_fun< RefEntity, EntityType, &RefEntity::getEntType >, const_mem_fun< RefEntity, EntityHandle, &RefEntity::getParentEnt > > > > > RefEntity_multiIndex
static auto get_sub_iface_options_imp(T *const ptr, int) -> decltype(ptr->getSubInterfaceOptions())
Definition Core.cpp:144
void set_ref_ent_basic_data_ptr_impl(boost::shared_ptr< BasicEntityData > &ptr)
Definition Core.cpp:921
multi_index_container< boost::shared_ptr< RefElement >, indexed_by< ordered_unique< tag< Ent_mi_tag >, const_mem_fun< RefElement::interface_type_RefEntity, EntityHandle, &RefElement::getEnt > > > > RefElement_multiIndex
static auto get_event_options_imp(T *const ptr, int) -> decltype(ptr->getEventOptions())
Definition Core.cpp:157
multi_index_container< boost::shared_ptr< FieldEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, member< FieldEntity, UId, &FieldEntity::localUId > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< FieldEntity::interface_type_RefEntity, EntityHandle, &FieldEntity::getEnt > > > > FieldEntity_multiIndex
static MoFEMErrorCode fixTagSize(moab::Interface &moab, bool *changed=nullptr)
Fix tag size when BITREFLEVEL_SIZE of core library is different than file BITREFLEVEL_SIZE.
Core (interface) class.
Definition Core.hpp:82
MoFEMErrorCode clearMap()
Cleaning database.
Definition Core.cpp:525
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Getting interface of core database.
Definition Core.cpp:43
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition Core.cpp:72
MoFEMErrorCode regSubInterface()
Register sub-interfaces in core interface.
Definition Core.cpp:172
static bool isGloballyInitialised
Core base globally initialized.
Definition Core.hpp:1056
const RefElement_multiIndex * get_ref_finite_elements() const
Get the ref finite elements object.
Definition Core.cpp:899
const FieldEntity_multiIndex * get_field_ents() const
Get the field ents object.
Definition Core.cpp:908
MoFEMErrorCode initialiseDatabaseFromMesh(int verb=DEFAULT_VERBOSITY)
Initialize database getting information on mesh.
Definition Core.cpp:300
MoFEMErrorCode registerSubInterfaces()
Register insterfaces.
Definition Core.cpp:485
const Field_multiIndex * get_fields() const
Get the fields object.
Definition Core.cpp:895
const Problem_multiIndex * get_problems() const
Get the problems object.
Definition Core.cpp:918
MoFEMErrorCode coreGenericConstructor(moab::Interface &moab, MPI_Comm comm, const int verbose)
Definition Core.cpp:203
const DofEntity_multiIndex * get_dofs() const
Get the dofs object.
Definition Core.cpp:911
MoFEMErrorCode getOptions(int verb=DEFAULT_VERBOSITY)
Get core options from command line.
Definition Core.cpp:777
static PetscBool isInitialized
petsc was initialised by other agent
Definition Core.hpp:1058
const RefEntity_multiIndex * get_ref_ents() const
Get the ref ents object.
Definition Core.cpp:896
const FieldEntityEntFiniteElementAdjacencyMap_multiIndex * get_ents_elements_adjacency() const
Get the dofs elements adjacency object.
Definition Core.cpp:891
const EntFiniteElement_multiIndex * get_ents_finite_elements() const
Get the ents finite elements object.
Definition Core.cpp:905
const FiniteElement_multiIndex * get_finite_elements() const
Get the finite elements object.
Definition Core.cpp:902
MoFEMErrorCode addPrismToDatabase(const EntityHandle prism, int verb=DEFAULT_VERBOSITY)
add prim element
Definition Core.cpp:544
MoFEMErrorCode regEvents()
Register petsc events.
Definition Core.cpp:190
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition Core.cpp:118
MoFEMErrorCode set_moab_interface(moab::Interface &new_moab, int verb=VERBOSE)
Set the moab interface object.
Definition Core.cpp:768
MoFEMErrorCode setMoabInterface(moab::Interface &new_moab, int verb=VERBOSE)
Definition Core.cpp:455
MoFEMErrorCode get_problem(const std::string &problem_name, const Problem **problem_ptr) const
Get problem database (data structure)
Definition Core.cpp:822
MoFEMErrorCode getTags(int verb=DEFAULT_VERBOSITY)
Get tag handles.
Definition Core.cpp:573
CoreTmp(moab::Interface &moab, MPI_Comm comm=PETSC_COMM_WORLD, const int verbose=VERBOSE)
Definition Core.cpp:241
static int mpiInitialised
mpi was initialised by other agent
Definition Core.hpp:1057
MeshsetsManager * get_meshsets_manager_ptr()
get MeshsetsManager pointer
Definition Core.cpp:870
MoFEMErrorCode rebuild_database(int verb=DEFAULT_VERBOSITY)
Clear database and initialize it once again.
Definition Core.cpp:758
MoFEMErrorCode clear_database(int verb=DEFAULT_VERBOSITY)
Clear database.
Definition Core.cpp:750
MoFEMErrorCode set_moab_interface(moab::Interface &new_moab, int verb)
CoreTmp(moab::Interface &moab, MPI_Comm comm=PETSC_COMM_WORLD, const int verbose=VERBOSE)
Deprecated interface functions.
Finite element definition.
static MoFEMErrorCode getOptions()
Get logger option.
static void createDefaultSinks(MPI_Comm comm)
Create default sinks.
static PetscErrorCode logPetscFPrintf(FILE *fd, const char format[], va_list Argp)
Use to handle PETSc output.
Interface for managing meshsets containing materials and boundary conditions.
keeps basic data about problem
keeps data about abstract PRISM finite element
MoFEMErrorCode getTags(int verb=-1)
get tags handlers used on meshsets
base class for all interface classes
MPI_Comm duplicatedComm
Definition Core.hpp:27
WrapMPIComm(MPI_Comm comm, bool petsc)
Definition Core.cpp:16
MPI_Comm comm
Definition Core.hpp:26