v0.15.0
Loading...
Searching...
No Matches
Namespaces | Macros
ProblemsCore.cpp File Reference

Managing complexities for problem. More...

#include <MoFEM.hpp>

Go to the source code of this file.

Namespaces

namespace  MoFEM
 implementation of Data Operators for Forces and Sources
 

Macros

#define ProblemCoreFunctionBegin
 
#define SET_BASIC_METHOD(METHOD, PROBLEM_PTR)
 

Detailed Description

Managing complexities for problem.

Definition in file ProblemsCore.cpp.

Macro Definition Documentation

◆ ProblemCoreFunctionBegin

#define ProblemCoreFunctionBegin
Value:
MOFEM_LOG_CHANNEL("WORLD"); \
MOFEM_LOG_CHANNEL("SYNC"); \
MOFEM_LOG_FUNCTION(); \
MOFEM_LOG_TAG("SYNC", "ProblemCore"); \
MOFEM_LOG_TAG("WORLD", "ProblemCore")
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...

Definition at line 8 of file ProblemsCore.cpp.

15 {
16
17bool Core::check_problem(const string name) {
18 Problem_multiIndex::index<Problem_mi_tag>::type::iterator pit;
19 pit = pRoblems.get<Problem_mi_tag>().find(name);
20 if (pit == pRoblems.get<Problem_mi_tag>().end()) {
21 return false;
22 }
23 return true;
24}
25
26MoFEMErrorCode Core::addProblem(const BitProblemId id, const std::string &name,
27 int verb) {
29
30 if (verb == -1)
31 verb = verbose;
32 EntityHandle meshset;
33 CHKERR get_moab().create_meshset(MESHSET_SET, meshset);
34 CHKERR get_moab().tag_set_data(th_ProblemId, &meshset, 1, &id);
35
36 // Add problem meshset to partion meshset. In case of no elements
37 // on processor part, when mesh file is read, finite element meshset is
38 // prevented from deletion by moab reader.
39 auto add_meshset_to_partition = [&](auto meshset) {
41 const void *tag_vals[] = {&rAnk};
42 ParallelComm *pcomm = ParallelComm::get_pcomm(
43 &get_moab(), get_basic_entity_data_ptr()->pcommID);
44 Tag part_tag = pcomm->part_tag();
45 Range tagged_sets;
46 CHKERR get_moab().get_entities_by_type_and_tag(0, MBENTITYSET, &part_tag,
47 tag_vals, 1, tagged_sets,
48 moab::Interface::UNION);
49 for (auto s : tagged_sets)
50 CHKERR get_moab().add_entities(s, &meshset, 1);
52 };
53 CHKERR add_meshset_to_partition(meshset);
54
55 void const *tag_data[] = {name.c_str()};
56 int tag_sizes[1];
57 tag_sizes[0] = name.size();
58 CHKERR get_moab().tag_set_by_ptr(th_ProblemName, &meshset, 1, tag_data,
59 tag_sizes);
60 // create entry
61 auto p = pRoblems.insert(Problem(moab, meshset));
62 if (!p.second) {
64 MOFEM_LOG_ATTRIBUTES("SELF", LogManager::BitScope);
65 MOFEM_LOG("SELF", Sev::error) << "Following problem can not be added:";
66 MOFEM_LOG("SELF", Sev::error)
67 << "Problem " << name << " id(" << id.to_ulong() << ") " << id
68 << " added meshset " << meshset;
69 MOFEM_LOG("SELF", Sev::error) << "List of problems already in databse:";
70 for (auto &p : pRoblems) {
71 MOFEM_LOG("SELF", Sev::error)
72 << p.getName() << " id(" << p.getId().to_ulong() << ") " << id
73 << " on meshset " << p.meshset;
74 }
75 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Problem not added");
76 }
77
78 MOFEM_LOG("WORLD", Sev::inform) << "Add problem " << name;
79
81}
82
83MoFEMErrorCode Core::add_problem(const std::string &name, enum MoFEMTypes bh,
84 int verb) {
86 if (verb == -1)
87 verb = verbose;
88 auto miit = pRoblems.get<Problem_mi_tag>().find(name);
89 if (miit == pRoblems.get<Problem_mi_tag>().end()) {
90
91 int p_shift = 0;
92 for (; pRoblems.get<BitProblemId_mi_tag>().find(BitProblemId().set(
93 p_shift)) != pRoblems.get<BitProblemId_mi_tag>().end();
94 ++p_shift) {
95 }
96
97 auto id = BitProblemId().set(p_shift);
98 CHKERR addProblem(id, name, verb);
99
100 } else if (bh == MF_EXCL) {
101 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_FOUND, "problem is in database %s",
102 name.c_str());
103 }
105}
106
107MoFEMErrorCode Core::delete_problem(const std::string name) {
109 auto p_miit = pRoblems.get<Problem_mi_tag>().find(name);
110 if (p_miit == pRoblems.get<Problem_mi_tag>().end()) {
111 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_FOUND, "no such problem like < %s >",
112 name.c_str());
113 }
114 const EntityHandle meshset = p_miit->meshset;
115 pRoblems.get<Problem_mi_tag>().erase(p_miit);
116 CHKERR get_moab().delete_entities(&meshset, 1);
118}
119
120BitProblemId Core::getBitProblemId(const std::string &name) const {
121 auto p_miit = pRoblems.get<Problem_mi_tag>().find(name);
122 if (p_miit == pRoblems.get<Problem_mi_tag>().end()) {
123 THROW_MESSAGE("no such problem like " + name + " >");
124 }
125 return p_miit->getId();
126}
127
128MoFEMErrorCode Core::list_problem() const {
130 typedef Problem_multiIndex::index<BitProblemId_mi_tag>::type ProblemById;
131 const ProblemById &set_id = pRoblems.get<BitProblemId_mi_tag>();
132 ProblemById::iterator miit = set_id.begin();
133 for (; miit != set_id.end(); miit++) {
134 std::ostringstream ss;
135 ss << *miit << std::endl;
136 PetscPrintf(mofemComm, "%s", ss.str().c_str());
137 }
139}
140
142Core::modify_problem_add_finite_element(const std::string name_problem,
143 const std::string &fe_name) {
145 typedef Problem_multiIndex::index<Problem_mi_tag>::type ProblemsByName;
146 ProblemsByName &set = pRoblems.get<Problem_mi_tag>();
147 ProblemsByName::iterator miit = set.find(name_problem);
148 if (miit == set.end()) {
149 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_FOUND, "this problem <%s> is not there",
150 name_problem.c_str());
151 }
152 BitFEId f_id = getBitFEId(fe_name);
153 bool success = set.modify(miit, ProblemFiniteElementChangeBitAdd(f_id));
154 if (!success)
155 SETERRQ(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
156 "modification unsuccessful");
158}
159
161Core::modify_problem_unset_finite_element(const std::string name_problem,
162 const std::string &fe_name) {
164 typedef Problem_multiIndex::index<Problem_mi_tag>::type ProblemsByName;
165 ProblemsByName &set = pRoblems.get<Problem_mi_tag>();
166 ProblemsByName::iterator miit = set.find(name_problem);
167 if (miit == set.end()) {
168 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_FOUND, "this problem <%s> is not there",
169 name_problem.c_str());
170 }
171 BitFEId f_id = getBitFEId(fe_name);
172 bool success = set.modify(miit, ProblemFiniteElementChangeBitUnSet(f_id));
173 if (!success)
174 SETERRQ(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
175 "modification unsuccessful");
177}
178
180Core::modify_problem_ref_level_add_bit(const std::string &name_problem,
181 const BitRefLevel &bit) {
183 typedef Problem_multiIndex::index<Problem_mi_tag>::type ProblemsByName;
184 ProblemsByName &set = pRoblems.get<Problem_mi_tag>();
185 ProblemsByName::iterator miit = set.find(name_problem);
186 std::ostringstream ss;
187 ss << name_problem;
188 if (miit == set.end())
189 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_FOUND, "this problem <%s> is there",
190 ss.str().c_str());
191 bool success = set.modify(miit, ProblemChangeRefLevelBitAdd(bit));
192 if (!success)
193 SETERRQ(PETSC_COMM_SELF, 1, "modification unsuccessful");
195}
196
198Core::modify_problem_ref_level_set_bit(const std::string &name_problem,
199 const BitRefLevel &bit) {
201 typedef Problem_multiIndex::index<Problem_mi_tag>::type ProblemsByName;
202 ProblemsByName &set = pRoblems.get<Problem_mi_tag>();
203 ProblemsByName::iterator miit = set.find(name_problem);
204 std::ostringstream ss;
205 ss << name_problem;
206 if (miit == set.end())
207 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_FOUND, "this problem <%s> is there",
208 ss.str().c_str());
209 bool success = set.modify(miit, ProblemChangeRefLevelBitSet(bit));
210 if (!success)
211 SETERRQ(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
212 "modification unsuccessful");
214}
215
217Core::modify_problem_mask_ref_level_add_bit(const std::string &name_problem,
218 const BitRefLevel &bit) {
220 typedef Problem_multiIndex::index<Problem_mi_tag>::type ProblemsByName;
221 ProblemsByName &set = pRoblems.get<Problem_mi_tag>();
222 ProblemsByName::iterator miit = set.find(name_problem);
223 if (miit == set.end()) {
224 std::ostringstream ss;
225 ss << name_problem;
226 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_FOUND, "this problem <%s> is there",
227 ss.str().c_str());
228 }
229 bool success = set.modify(miit, ProblemChangeRefLevelBitDofMaskAdd(bit));
230 if (!success)
231 SETERRQ(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
232 "modification unsuccessful");
234}
235
237Core::modify_problem_mask_ref_level_set_bit(const std::string &name_problem,
238 const BitRefLevel &bit) {
239 MoFEMFunctionBegin typedef Problem_multiIndex::index<Problem_mi_tag>::type
240 ProblemsByName;
241 ProblemsByName &set = pRoblems.get<Problem_mi_tag>();
242 ProblemsByName::iterator miit = set.find(name_problem);
243 if (miit == set.end()) {
244 std::ostringstream ss;
245 ss << name_problem;
246 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_FOUND, "this problem <%s> is there",
247 ss.str().c_str());
248 }
249 bool success = set.modify(miit, ProblemChangeRefLevelBitDofMaskSet(bit));
250 if (!success)
251 SETERRQ(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
252 "modification unsuccessful");
254}
255
256MoFEMErrorCode Core::build_problem_on_distributed_mesh(int verb) {
258 if (verb == -1)
259 verb = verbose;
260 Problem_multiIndex::iterator p_miit = pRoblems.begin();
261 for (; p_miit != pRoblems.end(); p_miit++) {
262 CHKERR getInterface<ProblemsManager>()->buildProblemOnDistributedMesh(
263 const_cast<Problem *>(&*p_miit), verb);
264 }
266}
267
268MoFEMErrorCode Core::clear_problem(const std::string problem_name, int verb) {
270 if (verb == -1)
271 verb = verbose;
272 typedef Problem_multiIndex::index<Problem_mi_tag>::type ProblemsByName;
273 ProblemsByName &prob_by_name = pRoblems.get<Problem_mi_tag>();
274 ProblemsByName::iterator p_miit = prob_by_name.find(problem_name);
275 if (p_miit == prob_by_name.end()) {
276 SETERRQ(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
277 "problem < %s > not found, (top tip: check spelling)",
278 problem_name.c_str());
279 }
280 // zero rows
281 bool success = prob_by_name.modify(p_miit, ProblemZeroNbRowsChange());
282 if (!success)
283 SETERRQ(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
284 "modification unsuccessful");
285 // zero cols
286 success = prob_by_name.modify(p_miit, ProblemZeroNbColsChange());
287 if (!success)
288 SETERRQ(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
289 "modification unsuccessful");
290 // clear finite elements
291 success =
292 prob_by_name.modify(p_miit, ProblemClearNumeredFiniteElementsChange());
293 if (!success)
294 SETERRQ(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
295 "modification unsuccessful");
296 // clear data structures
297 success = prob_by_name.modify(p_miit, ProblemClearSubProblemData());
298 if (!success)
299 SETERRQ(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
300 "modification unsuccessful");
301 success = prob_by_name.modify(p_miit, ProblemClearComposedProblemData());
302 if (!success)
303 SETERRQ(PETSC_COMM_SELF, MOFEM_OPERATION_UNSUCCESSFUL,
304 "modification unsuccessful");
305 if (p_miit->getRowDofsSequence())
306 p_miit->getRowDofsSequence()->clear();
307 if (p_miit->getColDofsSequence())
308 p_miit->getColDofsSequence()->clear();
309 if (p_miit->getSubData())
310 p_miit->getSubData().reset();
311 if (p_miit->getComposedProblemsData())
312 p_miit->getComposedProblemsData().reset();
314}
315
316MoFEMErrorCode Core::build_problems(int verb) {
318 if (verb == -1)
319 verb = verbose;
320 if (!((*buildMoFEM) & BUILD_FIELD))
321 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "fields not build");
322 if (!((*buildMoFEM) & BUILD_FE))
323 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "FEs not build");
324 if (!((*buildMoFEM) & BUILD_ADJ))
325 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "adjacencies not build");
326 // iterate problems
327 Problem_multiIndex::iterator p_miit = pRoblems.begin();
328 for (; p_miit != pRoblems.end(); p_miit++) {
329 Problem *problem_ptr = const_cast<Problem *>(&*p_miit);
330 CHKERR getInterface<ProblemsManager>()->buildProblem(problem_ptr, false,
331 verb);
332 }
333 *buildMoFEM |= BUILD_PROBLEM;
335}
336
337MoFEMErrorCode Core::clear_problems(int verb) {
339 if (verb == -1)
340 verb = verbose;
341 // iterate problems
342 for (auto p_miit = pRoblems.begin(); p_miit != pRoblems.end(); p_miit++)
343 CHKERR clear_problem(p_miit->getName(), verb);
345}
346
347#define SET_BASIC_METHOD(METHOD, PROBLEM_PTR) \
348 { \
349 METHOD.rAnk = rAnk; \
350 METHOD.sIze = sIze; \
351 METHOD.problemPtr = PROBLEM_PTR; \
352 METHOD.fieldsPtr = &fIelds; \
353 METHOD.refinedEntitiesPtr = &refinedEntities; \
354 METHOD.entitiesPtr = &entsFields; \
355 METHOD.dofsPtr = &dofsField; \
356 METHOD.refinedFiniteElementsPtr = &refinedFiniteElements; \
357 METHOD.finiteElementsPtr = &finiteElements; \
358 METHOD.finiteElementsEntitiesPtr = &entsFiniteElements; \
359 METHOD.adjacenciesPtr = &entFEAdjacencies; \
360 }
361
362MoFEMErrorCode Core::problem_basic_method_preProcess(const Problem *problem_ptr,
363 BasicMethod &method,
364 int verb) {
366 if (verb == -1)
367 verb = verbose;
368 // finite element
369 SET_BASIC_METHOD(method, problem_ptr)
370 PetscLogEventBegin(MOFEM_EVENT_preProcess, 0, 0, 0, 0);
371 CHKERR method.preProcess();
372 PetscLogEventEnd(MOFEM_EVENT_preProcess, 0, 0, 0, 0);
374}
375
377Core::problem_basic_method_preProcess(const std::string &problem_name,
378 BasicMethod &method, int verb) {
380 if (verb == -1)
381 verb = verbose;
382 typedef Problem_multiIndex::index<Problem_mi_tag>::type ProblemsByName;
383 // find p_miit
384 ProblemsByName &pRoblems_set = pRoblems.get<Problem_mi_tag>();
385 ProblemsByName::iterator p_miit = pRoblems_set.find(problem_name);
386 if (p_miit == pRoblems_set.end())
387 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_FOUND, "problem is not in database %s",
388 problem_name.c_str());
389 CHKERR problem_basic_method_preProcess(&*p_miit, method, verb);
391}
392
394Core::problem_basic_method_postProcess(const Problem *problem_ptr,
395 BasicMethod &method, int verb) {
397 SET_BASIC_METHOD(method, problem_ptr)
398
399 PetscLogEventBegin(MOFEM_EVENT_postProcess, 0, 0, 0, 0);
400 CHKERR method.postProcess();
401 PetscLogEventEnd(MOFEM_EVENT_postProcess, 0, 0, 0, 0);
402
404}
405
407Core::problem_basic_method_postProcess(const std::string &problem_name,
408 BasicMethod &method, int verb) {
410 if (verb == -1)
411 verb = verbose;
412 typedef Problem_multiIndex::index<Problem_mi_tag>::type ProblemsByName;
413
414 // find p_miit
415 ProblemsByName &pRoblems_set = pRoblems.get<Problem_mi_tag>();
416 ProblemsByName::iterator p_miit = pRoblems_set.find(problem_name);
417 if (p_miit == pRoblems_set.end())
418 SETERRQ(mofemComm, 1, "problem is not in database %s",
419 problem_name.c_str());
420
421 CHKERR problem_basic_method_postProcess(&*p_miit, method, verb);
422
424}
425
426MoFEMErrorCode Core::loop_finite_elements(
427 const Problem *problem_ptr, const std::string &fe_name, FEMethod &method,
428 int lower_rank, int upper_rank,
429 boost::shared_ptr<NumeredEntFiniteElement_multiIndex> fe_ptr, MoFEMTypes bh,
430 CacheTupleWeakPtr cache_ptr, int verb) {
432 if (verb == DEFAULT_VERBOSITY)
433 verb = verbose;
434
435 CacheTupleSharedPtr tmp_cache_ptr;
436 if (!cache_ptr.use_count()) {
437 tmp_cache_ptr = boost::make_shared<CacheTuple>();
438 CHKERR cache_problem_entities(problem_ptr->getName(), tmp_cache_ptr);
439 method.cacheWeakPtr = tmp_cache_ptr;
440 } else {
441 method.cacheWeakPtr = cache_ptr;
442 }
443
444 if (!fe_ptr)
445 fe_ptr = problem_ptr->numeredFiniteElementsPtr;
446
447 auto miit = fe_ptr->get<Composite_Name_And_Part_mi_tag>().lower_bound(
448 boost::make_tuple(fe_name, lower_rank));
449 auto hi_miit = fe_ptr->get<Composite_Name_And_Part_mi_tag>().upper_bound(
450 boost::make_tuple(fe_name, upper_rank));
451
452 if (miit == hi_miit && (bh & MF_EXIST)) {
453 if (!check_finite_element(fe_name)) {
454 SETERRQ(mofemComm, MOFEM_NOT_FOUND, "finite element < %s > not found",
455 fe_name.c_str());
456 }
457 }
458
459 method.feName = fe_name;
460 method.loopSize =
461 std::distance(miit, hi_miit); // Set numbers of elements in the loop
462 method.loHiFERank = std::make_pair(lower_rank, upper_rank);
463
464 SET_BASIC_METHOD(method, &*problem_ptr)
465
466 PetscLogEventBegin(MOFEM_EVENT_preProcess, 0, 0, 0, 0);
467 CHKERR method.preProcess();
468 PetscLogEventEnd(MOFEM_EVENT_preProcess, 0, 0, 0, 0);
469
470 PetscLogEventBegin(MOFEM_EVENT_operator, 0, 0, 0, 0);
471 for (int nn = 0; miit != hi_miit; miit++, nn++) {
472
473 method.nInTheLoop = nn; // Index of element in the loop
474 method.numeredEntFiniteElementPtr = *miit;
475
476 if (method.exeTestHook) {
477 if (method.exeTestHook(&method)) {
478 CHKERR method();
479 }
480 } else {
481 CHKERR method();
482 }
483
484 }
485 PetscLogEventEnd(MOFEM_EVENT_operator, 0, 0, 0, 0);
486
487 PetscLogEventBegin(MOFEM_EVENT_postProcess, 0, 0, 0, 0);
488 CHKERR method.postProcess();
489 PetscLogEventEnd(MOFEM_EVENT_postProcess, 0, 0, 0, 0);
490
492}
493
494MoFEMErrorCode Core::loop_finite_elements(
495 const std::string problem_name, const std::string &fe_name,
496 FEMethod &method,
497 boost::shared_ptr<NumeredEntFiniteElement_multiIndex> fe_ptr, MoFEMTypes bh,
498 CacheTupleWeakPtr cache_ptr, int verb) {
500 if (verb == DEFAULT_VERBOSITY)
501 verb = verbose;
502
503 CHKERR loop_finite_elements(problem_name, fe_name, method, rAnk, rAnk, fe_ptr,
504 bh, cache_ptr, verb);
505
507}
508
509MoFEMErrorCode Core::loop_finite_elements(
510 const std::string problem_name, const std::string &fe_name,
511 FEMethod &method, int lower_rank, int upper_rank,
512 boost::shared_ptr<NumeredEntFiniteElement_multiIndex> fe_ptr, MoFEMTypes bh,
513 CacheTupleWeakPtr cache_ptr, int verb) {
515 if (verb == DEFAULT_VERBOSITY)
516 verb = verbose;
517
518 auto &prb_by_name = pRoblems.get<Problem_mi_tag>();
519 auto p_miit = prb_by_name.find(problem_name);
520 if (p_miit == prb_by_name.end())
521 SETERRQ(mofemComm, MOFEM_INVALID_DATA, "Problem <%s> is not in database",
522 problem_name.c_str());
523
524 CHKERR loop_finite_elements(&*p_miit, fe_name, method, lower_rank, upper_rank,
525 fe_ptr, bh, cache_ptr, verb);
526
528}
529
530MoFEMErrorCode Core::loop_dofs(const Problem *problem_ptr,
531 const std::string &field_name, RowColData rc,
532 DofMethod &method, int lower_rank,
533 int upper_rank, int verb) {
535 SET_BASIC_METHOD(method, &*problem_ptr);
536 typedef NumeredDofEntity_multiIndex::index<Unique_mi_tag>::type
537 NumeredDofsByUId;
538 NumeredDofsByUId *dofs;
539 switch (rc) {
540 case ROW:
541 dofs = &problem_ptr->numeredRowDofsPtr->get<Unique_mi_tag>();
542 break;
543 case COL:
544 dofs = &problem_ptr->numeredColDofsPtr->get<Unique_mi_tag>();
545 break;
546 default:
547 SETERRQ(mofemComm, MOFEM_DATA_INCONSISTENCY, "Not implemented");
548 }
549
550 auto field_it = fIelds.get<FieldName_mi_tag>().find(field_name);
551 if (field_it != fIelds.get<FieldName_mi_tag>().end()) {
552 method.fieldPtr = *field_it;
553 } else {
554 SETERRQ(mofemComm, MOFEM_NOT_FOUND, "Field not found %s",
555 field_name.c_str());
556 }
557
558 auto miit = dofs->lower_bound(
559 FieldEntity::getLoBitNumberUId((*field_it)->getBitNumber()));
560 auto hi_miit = dofs->upper_bound(
561 FieldEntity::getHiBitNumberUId((*field_it)->getBitNumber()));
562
563 method.loopSize = std::distance(miit, hi_miit);
564 method.loHiFERank = std::make_pair(lower_rank, upper_rank);
565
566 CHKERR method.preProcess();
567
568 int nn = 0;
569 for (; miit != hi_miit; miit++, nn++) {
570 if ((*miit)->getPart() >= lower_rank && (*miit)->getPart() <= upper_rank) {
571 method.nInTheLoop = nn; // Index of element in the loop
572 method.dofPtr = miit->get()->getDofEntityPtr();
573 method.dofNumeredPtr = *miit;
574 CHKERR method();
575 }
576 }
577
578 CHKERR method.postProcess();
580}
581
582MoFEMErrorCode Core::loop_dofs(
583 const std::string &problem_name, const std::string &field_name,
584 RowColData rc, // ROW or COL
585 DofMethod &method, // Finite element instance processed on each DOF
586 int lower_rank, // Only DOFs on processor higher or equal to this are
587 // processed
588 int upper_rank, // Only DOFs lowest or higher to this are processed
589 int verb // verbosity level
590) {
592 if (verb == DEFAULT_VERBOSITY)
593 verb = verbose;
594 typedef Problem_multiIndex::index<Problem_mi_tag>::type ProblemsByName;
595 // find p_miit
596 ProblemsByName &pRoblems_set = pRoblems.get<Problem_mi_tag>();
597 ProblemsByName::iterator p_miit = pRoblems_set.find(problem_name);
598 if (p_miit == pRoblems_set.end())
599 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_FOUND, "problem not in database %s",
600 problem_name.c_str());
601 CHKERR loop_dofs(&*p_miit, field_name, rc, method, lower_rank, upper_rank,
602 verb);
604}
605
606MoFEMErrorCode Core::loop_dofs(const std::string &problem_name,
607 const std::string &field_name, RowColData rc,
608 DofMethod &method, int verb) {
610 if (verb == DEFAULT_VERBOSITY)
611 verb = verbose;
612 CHKERR loop_dofs(problem_name, field_name, rc, method, 0, sIze, verb);
614}
615
616MoFEMErrorCode Core::loop_dofs(const std::string &field_name, DofMethod &method,
617 int verb) {
619 if (verb == DEFAULT_VERBOSITY)
620 verb = verbose;
621 SET_BASIC_METHOD(method, nullptr);
622
623 auto field_it = fIelds.get<FieldName_mi_tag>().find(field_name);
624 if (field_it != fIelds.get<FieldName_mi_tag>().end()) {
625 method.fieldPtr = *field_it;
626 } else {
627 SETERRQ(mofemComm, MOFEM_NOT_FOUND, "Field not found %s",
628 field_name.c_str());
629 }
630
631 auto miit = dofsField.get<Unique_mi_tag>().lower_bound(
632 FieldEntity::getLoBitNumberUId((*field_it)->getBitNumber()));
633 auto hi_miit = dofsField.get<Unique_mi_tag>().upper_bound(
634 FieldEntity::getHiBitNumberUId((*field_it)->getBitNumber()));
635
636 method.loopSize = std::distance(miit, hi_miit);
637 method.loHiFERank = std::make_pair(0, get_comm_size());
638
639 CHKERR method.preProcess();
640 for (int nn = 0; miit != hi_miit; miit++, nn++) {
641 method.nInTheLoop = nn;
642 method.dofPtr = *miit;
643 CHKERR method();
644 }
645 CHKERR method.postProcess();
647}
648
649MoFEMErrorCode Core::loop_entities(const Problem *problem_ptr,
650 const std::string field_name, RowColData rc,
651 EntityMethod &method, int lower_rank,
652 int upper_rank, int verb) {
654 if (verb == DEFAULT_VERBOSITY)
655 verb = verbose;
656 decltype(problem_ptr->numeredRowDofsPtr) dofs;
657 switch (rc) {
658 case ROW:
659 dofs = problem_ptr->numeredRowDofsPtr;
660 break;
661 case COL:
662 dofs = problem_ptr->numeredColDofsPtr;
663 break;
664 default:
665 SETERRQ(mofemComm, MOFEM_DATA_INCONSISTENCY,
666 "It works only with rows or columns");
667 }
668
669 auto field_it = fIelds.get<FieldName_mi_tag>().find(field_name);
670 if (field_it != fIelds.get<FieldName_mi_tag>().end()) {
671 method.fieldPtr = *field_it;
672 } else {
673 SETERRQ(mofemComm, MOFEM_NOT_FOUND, "Field not found %s",
674 field_name.c_str());
675 }
676
677 auto miit = dofs->lower_bound(
678 FieldEntity::getLoBitNumberUId((*field_it)->getBitNumber()));
679 auto hi_miit = dofs->upper_bound(
680 FieldEntity::getHiBitNumberUId((*field_it)->getBitNumber()));
681
682 using FieldEntity_view_multiIndex = multi_index_container<
683
684 boost::shared_ptr<FieldEntity>,
685 indexed_by<
686
687 ordered_unique<
688
689 tag<Ent_mi_tag>,
690 const_mem_fun<FieldEntity::interface_type_RefEntity, EntityHandle,
691 &FieldEntity::getEnt>>
692
693 >>;
694
695 FieldEntity_view_multiIndex ents_view;
696 auto hint = ents_view.begin();
697 for (; miit != hi_miit; ++miit)
698 if ((*miit)->getPart() >= lower_rank && (*miit)->getPart() <= upper_rank)
699 ents_view.emplace_hint(hint, (*miit)->getFieldEntityPtr());
700
701 SET_BASIC_METHOD(method, problem_ptr);
702
703 method.loopSize = ents_view.size();
704 method.loHiFERank = std::make_pair(lower_rank, upper_rank);
705
706 CHKERR method.preProcess();
707 method.nInTheLoop = 0;
708 for (auto &field_ent : ents_view) {
709 method.entPtr = field_ent;
710 CHKERR method();
711 ++method.nInTheLoop;
712 }
713 CHKERR method.postProcess();
715}
716
717MoFEMErrorCode Core::loop_entities(const std::string problem_name,
718 const std::string field_name, RowColData rc,
719 EntityMethod &method, int lower_rank,
720 int upper_rank, int verb) {
722 if (verb == DEFAULT_VERBOSITY)
723 verb = verbose;
724 // find p_miit
725 auto &prb = pRoblems.get<Problem_mi_tag>();
726 auto p_miit = prb.find(problem_name);
727 if (p_miit == prb.end())
728 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_FOUND, "problem not in database %s",
729 problem_name.c_str());
730 CHKERR loop_entities(&*p_miit, field_name, rc, method, lower_rank, upper_rank,
731 verb);
733}
734
735MoFEMErrorCode Core::loop_entities(const std::string problem_name,
736 const std::string field_name, RowColData rc,
737 EntityMethod &method, int verb) {
738 return loop_entities(problem_name, field_name, rc, method, rAnk, rAnk, verb);
739}
740
741MoFEMErrorCode Core::loop_entities(const std::string field_name,
742 EntityMethod &method,
743 Range const *const ents, int verb) {
745 if (verb == DEFAULT_VERBOSITY)
746 verb = verbose;
747 SET_BASIC_METHOD(method, nullptr);
748
749 auto field_it = fIelds.get<FieldName_mi_tag>().find(field_name);
750 if (field_it != fIelds.get<FieldName_mi_tag>().end()) {
751 method.fieldPtr = *field_it;
752 } else {
753 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_FOUND, "Field not found %s",
754 field_name.c_str());
755 }
756
757 auto lo_eit = entsFields.get<Unique_mi_tag>().lower_bound(
758 FieldEntity::getLoBitNumberUId((*field_it)->getBitNumber()));
759 auto hi_eit = entsFields.get<Unique_mi_tag>().upper_bound(
760 FieldEntity::getHiBitNumberUId((*field_it)->getBitNumber()));
761
762 typedef multi_index_container<
763 boost::shared_ptr<FieldEntity>,
764 indexed_by<ordered_unique<
765 tag<Ent_mi_tag>, const_mem_fun<FieldEntity::interface_RefEntity,
766 EntityHandle, &FieldEntity::getEnt>>>>
767 FieldEntity_view_multiIndex;
768
769 FieldEntity_view_multiIndex ents_view;
770 ents_view.insert(lo_eit, hi_eit);
771
772 method.loopSize = ents_view.size();
773 method.loHiFERank = std::make_pair(0, get_comm_size());
774
775 CHKERR method.preProcess();
776 method.nInTheLoop = 0;
777
778 if (ents)
779 for (auto p = ents->const_pair_begin(); p != ents->const_pair_end(); ++p)
780 for (auto feit = ents_view.lower_bound(p->first);
781 feit != ents_view.upper_bound(p->second); ++feit) {
782 method.entPtr = *feit;
783 CHKERR method();
784 ++method.nInTheLoop;
785 }
786 else
787 for (auto &field_ent : ents_view) {
788 method.entPtr = field_ent;
789 CHKERR method();
790 ++method.nInTheLoop;
791 }
792
793 CHKERR method.postProcess();
795}
796
797MoFEMErrorCode Core::cache_problem_entities(const std::string prb_name,
798 CacheTupleWeakPtr cache_weak_ptr) {
800
801 if (auto cache_ptr = cache_weak_ptr.lock()) {
802 auto p_miit = pRoblems.get<Problem_mi_tag>().find(prb_name);
803 if (p_miit == pRoblems.get<Problem_mi_tag>().end())
804 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_FOUND, "problem not in database %s",
805 prb_name.c_str());
806
807 const BitRefLevel &prb_bit = p_miit->getBitRefLevel();
808 const BitRefLevel &prb_mask = p_miit->getBitRefLevelMask();
809 const BitFEId &prb_fe_id = p_miit->getBitFEId();
810
811 auto &row_dofs = p_miit->numeredRowDofsPtr;
812 auto &col_dofs = p_miit->numeredColDofsPtr;
813
814 auto &cache_data = std::get<0>(*cache_ptr);
815 auto &cache_row = std::get<1>(*cache_ptr);
816 auto &cache_col = std::get<2>(*cache_ptr);
817
818 cache_row.resize(entsFields.size());
819 if (row_dofs != col_dofs)
820 cache_col.resize(entsFields.size());
821 cache_data.resize(entsFields.size());
822
823 size_t idx = 0;
824 for (auto it = entsFields.begin(); it != entsFields.end(); ++it, ++idx) {
825
826 const auto uid = (*it)->getLocalUniqueId();
827 auto r = entFEAdjacencies.get<Unique_mi_tag>().equal_range(uid);
828 for (auto lo = r.first; lo != r.second; ++lo) {
829
830 if ((lo->getBitFEId() & prb_fe_id).any()) {
831
832 const BitRefLevel &fe_bit = lo->entFePtr->getBitRefLevel();
833
834 // if entity is not problem refinement level
835 if (((fe_bit & prb_mask) != fe_bit) || ((fe_bit & prb_bit).none()))
836 continue;
837
838 auto cache_numered_dofs = [&](auto &numered_dofs, auto &cache_vec,
839 auto &ent_cache) {
840 auto dit = numered_dofs->lower_bound(uid);
841 decltype(dit) hi_dit;
842 if (dit != numered_dofs->end())
843 hi_dit = numered_dofs->upper_bound(
844 uid | static_cast<UId>(MAX_DOFS_ON_ENTITY - 1));
845 else
846 hi_dit = dit;
847
848 ent_cache = boost::shared_ptr<EntityCacheNumeredDofs>(
849 cache_ptr, &(cache_vec[idx]));
850 cache_vec[idx].loHi = {dit, hi_dit};
851 };
852
853 auto cache_dofs = [&](auto &dofs, auto &cache_vec, auto &ent_cache) {
854 auto dit = dofs.lower_bound(uid);
855 decltype(dit) hi_dit;
856 if (dit != dofs.end())
857 hi_dit = dofs.upper_bound(
858 uid | static_cast<UId>(MAX_DOFS_ON_ENTITY - 1));
859 else
860 hi_dit = dit;
861
862 ent_cache = boost::shared_ptr<EntityCacheDofs>(cache_ptr,
863 &(cache_vec[idx]));
864 cache_vec[idx].loHi = {dit, hi_dit};
865 };
866
867 cache_numered_dofs(row_dofs, cache_row, (*it)->entityCacheRowDofs);
868 if (row_dofs != col_dofs) {
869 if (cache_col.size() != entsFields.size())
870 cache_col.resize(entsFields.size());
871 cache_numered_dofs(col_dofs, cache_col, (*it)->entityCacheColDofs);
872 } else {
873 (*it)->entityCacheColDofs = (*it)->entityCacheRowDofs;
874 }
875
876 cache_dofs(dofsField, cache_data, (*it)->entityCacheDataDofs);
877
878 break;
879 }
880 }
881 }
882 } else {
883 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Cache not allocated");
884 }
885
887}
888
889} // namespace MoFEM
#define ProblemCoreFunctionBegin
#define SET_BASIC_METHOD(METHOD, PROBLEM_PTR)
@ DEFAULT_VERBOSITY
RowColData
RowColData.
@ COL
@ ROW
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
@ MF_EXIST
@ MF_EXCL
#define MAX_DOFS_ON_ENTITY
Maximal number of DOFs on entity.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
@ MOFEM_NOT_FOUND
Definition definitions.h:33
@ MOFEM_OPERATION_UNSUCCESSFUL
Definition definitions.h:34
@ MOFEM_DATA_INCONSISTENCY
Definition definitions.h:31
@ MOFEM_INVALID_DATA
Definition definitions.h:36
#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.
#define MOFEM_LOG(channel, severity)
Log.
#define MOFEM_LOG_ATTRIBUTES(channel, bit)
Add attributes to channel.
#define MOFEM_LOG_FUNCTION()
Set scope.
auto bit
set bit
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
std::bitset< BITFEID_SIZE > BitFEId
Finite element Id.
Definition Types.hpp:43
std::bitset< BITPROBLEMID_SIZE > BitProblemId
Problem Id.
Definition Types.hpp:44
uint128_t UId
Unique Id.
Definition Types.hpp:31
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition Types.hpp:40
Field_multiIndex::index< FieldName_mi_tag >::type::iterator field_it
boost::shared_ptr< CacheTuple > CacheTupleSharedPtr
int r
Definition sdf.py:8
constexpr auto field_name

◆ SET_BASIC_METHOD

#define SET_BASIC_METHOD (   METHOD,
  PROBLEM_PTR 
)
Value:
{ \
METHOD.rAnk = rAnk; \
METHOD.sIze = sIze; \
METHOD.problemPtr = PROBLEM_PTR; \
METHOD.fieldsPtr = &fIelds; \
METHOD.refinedEntitiesPtr = &refinedEntities; \
METHOD.entitiesPtr = &entsFields; \
METHOD.dofsPtr = &dofsField; \
METHOD.refinedFiniteElementsPtr = &refinedFiniteElements; \
METHOD.finiteElementsPtr = &finiteElements; \
METHOD.finiteElementsEntitiesPtr = &entsFiniteElements; \
METHOD.adjacenciesPtr = &entFEAdjacencies; \
}

Definition at line 348 of file ProblemsCore.cpp.

349 { \
350 METHOD.rAnk = rAnk; \
351 METHOD.sIze = sIze; \
352 METHOD.problemPtr = PROBLEM_PTR; \
353 METHOD.fieldsPtr = &fIelds; \
354 METHOD.refinedEntitiesPtr = &refinedEntities; \
355 METHOD.entitiesPtr = &entsFields; \
356 METHOD.dofsPtr = &dofsField; \
357 METHOD.refinedFiniteElementsPtr = &refinedFiniteElements; \
358 METHOD.finiteElementsPtr = &finiteElements; \
359 METHOD.finiteElementsEntitiesPtr = &entsFiniteElements; \
360 METHOD.adjacenciesPtr = &entFEAdjacencies; \
361 }