v0.13.1
BitRefManager.cpp
Go to the documentation of this file.
1/** \file BitRefManager.cpp
2 * \brief Managing BitRefLevels
3 * \mofem_bit_ref
4 */
5
6/* MoFEM is free software: you can redistribute it and/or modify it under
7 * the terms of the GNU Lesser General Public License as published by the
8 * Free Software Foundation, either version 3 of the License, or (at your
9 * option) any later version.
10 *
11 * MoFEM is distributed in the hope that it will be useful, but WITHOUT
12 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
14 * License for more details.
15 *
16 * You should have received a copy of the GNU Lesser General Public
17 * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>
18 */
19
20namespace MoFEM {
21
23BitRefManager::query_interface(boost::typeindex::type_index type_index,
24 UnknownInterface **iface) const {
26 *iface = const_cast<BitRefManager *>(this);
28}
29
31 : cOre(const_cast<MoFEM::Core &>(core)), dEbug(false) {
32
33 if (!LogManager::checkIfChannelExist("BitRefSelf")) {
34 auto core_log = logging::core::get();
35 core_log->add_sink(
37 LogManager::setLog("BitRefSelf");
38 core_log->add_sink(
40 LogManager::setLog("BitRefWorld");
41 core_log->add_sink(
43 LogManager::setLog("BitRefSync");
44 MOFEM_LOG_TAG("BitRefSelf", "BitRefManager");
45 MOFEM_LOG_TAG("BitRefWorld", "BitRefManager");
46 MOFEM_LOG_TAG("BitRefSync", "BitRefManager");
47 }
48
50 MOFEM_LOG("BitRefWorld", Sev::noisy) << "BitRefManager interface created";
51}
52
53/// tool class with methods used more than twp times
55
57
58 const BitRefLevel &bIt; ///< bit to set
59 const RefEntity_multiIndex *refEntsPtr; ///< access to ents database
60 const RefElement_multiIndex *refElementPtr; ///< access to fe database
61
62 boost::shared_ptr<BasicEntityData> &baseEntData; ///< base entity data
63
64 /// constrictor
66 const RefEntity_multiIndex *ref_ents_ptr,
67 const RefElement_multiIndex *ref_element_ptr)
68 : mField(m_field), bIt(bit), refEntsPtr(ref_ents_ptr),
69 refElementPtr(ref_element_ptr),
70 baseEntData(m_field.get_basic_entity_data_ptr()) {}
71
72 /// find entities and change entity bit if in database
73 MoFEMErrorCode findEntsToAdd(EntityHandle f, EntityHandle s,
74 Range &seed_ents_range) const {
76
77 seed_ents_range.insert(f, s);
78 // get lower bound of multi-index
79 auto rit = refEntsPtr->lower_bound(f);
80 if (rit == refEntsPtr->end()) {
81 // all enties in range are added to database
83 } else {
84
85 // some entities from range are in database
86 auto hi_rit = refEntsPtr->upper_bound(s);
87
88 Range to_erase;
89 insertOrdered(to_erase, RefEntExtractor(), rit, hi_rit);
90 if (bIt.any())
91 for (; rit != hi_rit; ++rit)
92 const_cast<BitRefLevel &>((*rit)->getBitRefLevel()) |= bIt;
93
94 Range::iterator lo, hi = seed_ents_range.begin();
95 for (auto pt = to_erase.pair_begin(); pt != to_erase.pair_end(); ++pt) {
96 lo = seed_ents_range.lower_bound(hi, seed_ents_range.end(), pt->first);
97 if (lo != seed_ents_range.end()) {
98 hi = seed_ents_range.upper_bound(lo, seed_ents_range.end(),
99 pt->second);
100 seed_ents_range.erase(lo, hi);
101 } else
102 break;
103 }
104 }
105
107 }
108
109 /// add entities to database
110 template <int N>
111 MoFEMErrorCode addEntsToDatabaseImpl(const Range &seed_ents_range) const {
113 std::vector<boost::shared_ptr<RefEntity>> shared_ref_ents_vec;
114 shared_ref_ents_vec.reserve(seed_ents_range.size());
115 std::vector<const void *> tag_by_ptr;
116 for (Range::const_pair_iterator pit = seed_ents_range.pair_begin();
117 pit != seed_ents_range.pair_end(); pit++) {
118 // add entities to database
119 EntityHandle f = pit->first;
120 EntityHandle s = pit->second;
121
122 boost::shared_ptr<std::vector<RefEntityTmp<N>>> ref_ents_vec(
123 new std::vector<RefEntityTmp<N>>());
124 ref_ents_vec->reserve(s - f + 1);
125
126 tag_by_ptr.resize(s - f + 1);
127 CHKERR baseEntData->moab.tag_get_by_ptr(
128 baseEntData->th_RefParentHandle, Range(f, s), &*tag_by_ptr.begin());
129 auto tag_parent_it = tag_by_ptr.begin();
130 for (auto f : Range(f, s)) {
131 ref_ents_vec->emplace_back(
132 baseEntData, f,
133 const_cast<EntityHandle *>(
134 static_cast<const EntityHandle *>(*tag_parent_it)));
135 ++tag_parent_it;
136 }
137
138 // Set bits to range
139 if (bIt.any()) {
140 tag_by_ptr.resize(s - f + 1);
141 CHKERR baseEntData->moab.tag_get_by_ptr(
142 baseEntData->th_RefBitLevel, Range(f, s), &*tag_by_ptr.begin());
143 for (auto &v_bit_ptr : tag_by_ptr)
144 const_cast<BitRefLevel &>(
145 *(static_cast<const BitRefLevel *>(v_bit_ptr))) |= bIt;
146 }
147
148 for (auto &re : *ref_ents_vec)
149 shared_ref_ents_vec.emplace_back(ref_ents_vec,
150 static_cast<RefEntity *>(&re));
151 }
152 if (!shared_ref_ents_vec.empty()) {
153 int s0 = refEntsPtr->size();
154 const_cast<RefEntity_multiIndex *>(refEntsPtr)
155 ->insert(shared_ref_ents_vec.begin(), shared_ref_ents_vec.end());
156 if ((refEntsPtr->size() - s0) != shared_ref_ents_vec.size()) {
157 SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
158 "Data inconsistency %d != %d", refEntsPtr->size() - s0,
159 shared_ref_ents_vec.size());
160 }
161 }
163 }
164
165 MoFEMErrorCode addEntsToDatabase(const Range &seed_ents_range) const {
167
168 switch (mField.getValue()) {
169 case -1:
170 return addEntsToDatabaseImpl<-1>(seed_ents_range);
171 case 0:
172 return addEntsToDatabaseImpl<0>(seed_ents_range);
173 case 1:
174 return addEntsToDatabaseImpl<1>(seed_ents_range);
175 default:
176 SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
177 "Core index can vary from -1 to %d", MAX_CORE_TMP);
178 }
179
181 }
182
183 MoFEMErrorCode findElementsToAdd(EntityHandle f, EntityHandle s,
184 Range &seed_fe_range) const {
186 seed_fe_range.insert(f, s);
187 RefElement_multiIndex::iterator rit, hi_rit;
188 // get lower bound of multi-index
189 rit = refElementPtr->lower_bound(f);
190 if (rit == refElementPtr->end()) {
191 // all enties in range are added to database
193 } else {
194 // some entities from range are in database
195 hi_rit = refElementPtr->upper_bound(s);
196 for (; rit != hi_rit; ++rit) {
197 seed_fe_range.erase(rit->get()->getEnt());
198 }
199 }
201 }
202
203 MoFEMErrorCode addElementsToDatabase(Range &seed_fe_range) const {
205 std::vector<boost::shared_ptr<RefElement>> shared_ref_fe_vec;
206 shared_ref_fe_vec.reserve(seed_fe_range.size());
207 // create ref entity instances
208 for (Range::const_pair_iterator pit = seed_fe_range.const_pair_begin();
209 pit != seed_fe_range.const_pair_end(); ++pit) {
210 RefEntity_multiIndex::iterator rit, hi_rit;
211 rit = refEntsPtr->lower_bound(pit->first);
212 hi_rit = refEntsPtr->upper_bound(pit->second);
213 if (static_cast<int>(std::distance(rit, hi_rit)) !=
214 static_cast<int>(pit->second - pit->first + 1)) {
215 SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
216 "data inconsistency %d != %d", std::distance(rit, hi_rit),
217 pit->second - pit->first + 1);
218 }
219 switch ((*rit)->getEntType()) {
220 case MBVERTEX: {
221 boost::shared_ptr<std::vector<RefElement_VERTEX>> ref_fe_vec =
222 boost::make_shared<std::vector<RefElement_VERTEX>>();
223 ref_fe_vec->reserve(pit->second - pit->first + 1);
224 for (; rit != hi_rit; ++rit) {
225 ref_fe_vec->push_back(RefElement_VERTEX(*rit));
226 shared_ref_fe_vec.push_back(
227 boost::shared_ptr<RefElement>(ref_fe_vec, &ref_fe_vec->back()));
228 }
229 } break;
230 case MBEDGE: {
231 boost::shared_ptr<std::vector<RefElement_EDGE>> ref_fe_vec =
232 boost::make_shared<std::vector<RefElement_EDGE>>();
233 ref_fe_vec->reserve(pit->second - pit->first + 1);
234 for (; rit != hi_rit; ++rit) {
235 ref_fe_vec->push_back(RefElement_EDGE(*rit));
236 shared_ref_fe_vec.push_back(
237 boost::shared_ptr<RefElement>(ref_fe_vec, &ref_fe_vec->back()));
238 }
239 } break;
240 case MBTRI:
241 case MBQUAD: {
242 boost::shared_ptr<std::vector<RefElementFace>> ref_fe_vec =
243 boost::make_shared<std::vector<RefElementFace>>();
244 ref_fe_vec->reserve(pit->second - pit->first + 1);
245 for (; rit != hi_rit; ++rit) {
246 ref_fe_vec->push_back(RefElementFace(*rit));
247 shared_ref_fe_vec.push_back(
248 boost::shared_ptr<RefElement>(ref_fe_vec, &ref_fe_vec->back()));
249 }
250 } break;
251 case MBTET:
252 case MBHEX: {
253 boost::shared_ptr<std::vector<RefElementVolume>> ref_fe_vec =
254 boost::make_shared<std::vector<RefElementVolume>>();
255 ref_fe_vec->reserve(pit->second - pit->first + 1);
256 for (; rit != hi_rit; ++rit) {
257 ref_fe_vec->push_back(RefElementVolume(*rit));
258 shared_ref_fe_vec.push_back(
259 boost::shared_ptr<RefElement>(ref_fe_vec, &ref_fe_vec->back()));
260 }
261 } break;
262 case MBPRISM: {
263 boost::shared_ptr<std::vector<RefElement_PRISM>> ref_fe_vec =
264 boost::make_shared<std::vector<RefElement_PRISM>>();
265 ref_fe_vec->reserve(pit->second - pit->first + 1);
266 for (; rit != hi_rit; ++rit) {
267 ref_fe_vec->push_back(RefElement_PRISM(*rit));
268 shared_ref_fe_vec.push_back(
269 boost::shared_ptr<RefElement>(ref_fe_vec, &ref_fe_vec->back()));
270 }
271 } break;
272 case MBENTITYSET: {
273 boost::shared_ptr<std::vector<RefElement_MESHSET>> ref_fe_vec =
274 boost::make_shared<std::vector<RefElement_MESHSET>>();
275 ref_fe_vec->reserve(pit->second - pit->first + 1);
276 for (; rit != hi_rit; ++rit) {
277 ref_fe_vec->push_back(RefElement_MESHSET(*rit));
278 shared_ref_fe_vec.push_back(
279 boost::shared_ptr<RefElement>(ref_fe_vec, &ref_fe_vec->back()));
280 }
281 } break;
282 default:
283 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
284 }
285 }
286 // add shared pointers to database
288 ->insert(shared_ref_fe_vec.begin(), shared_ref_fe_vec.end());
290 }
291};
292
294 const BitRefLevel bit,
295 const bool only_tets, int verb,
296 Range *adj_ents_ptr) const {
297 MoFEM::Interface &m_field = cOre;
298 auto ref_ents_ptr = m_field.get_ref_ents();
299 auto ref_fe_ptr = m_field.get_ref_finite_elements();
301
303 MOFEM_LOG_C("BitRefSelf", Sev::noisy, "Number of entities to add %d",
304 ents.size());
305
306 CHKERR setElementsBitRefLevel(ents, bit, verb);
307
308 if (!ents.empty()) {
309 for (int d = 3; d >= 1; --d) {
310 Range dim_ents;
311 if (only_tets && d == 3) {
312 dim_ents = ents.subset_by_type(MBTET);
313 } else {
314 dim_ents = ents.subset_by_dimension(d);
315 }
316
318 MOFEM_LOG_C("BitRefSelf", Sev::noisy,
319 "\tNumber of dim %d entities to add %d", d, dim_ents.size());
320
321 if (!dim_ents.empty()) {
322 for (int dd = 0; dd < d; ++dd) {
323 Range adj_ents;
324
325 if (dd == 0) {
326 rval = m_field.get_moab().get_connectivity(ents, adj_ents, true);
327 // rval = m_field.get_moab().get_adjacencies(
328 // dim_ents, dd, true, adj_ents, moab::Interface::UNION);
329
330 } else {
331 if (adj_ents_ptr) {
332 if (dd == 1) {
333 adj_ents = adj_ents_ptr->subset_by_dimension(MBEDGE);
334 } else if (dd == 2) {
335 adj_ents = adj_ents_ptr->subset_by_dimension(MBTRI);
336 }
337 } else {
338 rval = m_field.get_moab().get_adjacencies(
339 dim_ents, dd, true, adj_ents, moab::Interface::UNION);
340 }
341 }
342
343 // rval = m_field.get_moab().get_adjacencies(
344 // dim_ents, dd, true, adj_ents, moab::Interface::UNION);
345
347 MOFEM_LOG_C("BitRefSelf", Sev::noisy,
348 "\tNumber of dim %d adj entities for dim %d to add %d", d,
349 dd, adj_ents.size());
350
351 if (rval == MB_MULTIPLE_ENTITIES_FOUND) {
352 auto log_message = [&](const auto sev) {
355 MOFEM_LOG("BitRefSelf", sev)
356 << "When get adjacencies moab return MB_MULTIPLE_ENTITIES_ "
357 "FOUND for dim = "
358 << dd << " and dim of entities " << d;
359 MOFEM_LOG_CHANNEL("BitRefSelf"); // reset channel
360 };
361
362 if (verb <= QUIET)
363 log_message(Sev::noisy);
364 else
365 log_message(Sev::warning);
366
367 rval = MB_SUCCESS;
368 }
370 for (Range::pair_iterator pit = adj_ents.pair_begin();
371 pit != adj_ents.pair_end(); ++pit) {
372 Range seed_ents_range;
373 // get first and last element of range
374 EntityHandle f = pit->first;
375 EntityHandle s = pit->second;
376 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
377 .findEntsToAdd(f, s, seed_ents_range);
378 if (!seed_ents_range.empty())
379 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
380 .addEntsToDatabase(seed_ents_range);
381 }
382 }
383 }
384 }
385 }
386
388}
389
391 const BitRefLevel bit,
392 int verb) const {
393 MoFEM::Interface &m_field = cOre;
394 auto ref_ents_ptr = m_field.get_ref_ents();
395 auto ref_fe_ptr = m_field.get_ref_finite_elements();
397
398 for (Range::const_pair_iterator pit = ents.pair_begin();
399 pit != ents.pair_end(); pit++) {
400 // get first and last element of range
401 EntityHandle f = pit->first;
402 EntityHandle s = pit->second;
403 Range seed_ents_range; // entities seeded not in database
404 // find ents to add
405 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
406 .findEntsToAdd(f, s, seed_ents_range);
407 // add elements
408 if (!seed_ents_range.empty())
409 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
410 .addEntsToDatabase(seed_ents_range);
411
412 Range seed_fe_range;
413 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
414 .findElementsToAdd(f, s, seed_fe_range);
415 if (!seed_fe_range.empty()) {
416 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
417 .addElementsToDatabase(seed_fe_range);
418 }
419 }
420
422 MOFEM_LOG("BitRefSelf", Sev::noisy)
423 << "Number of entities in databse " << ref_ents_ptr->size();
424 MOFEM_LOG("BitRefSelf", Sev::noisy)
425 << "Number of finite element entities in databse " << ref_fe_ptr->size();
426
428}
429
431 const BitRefLevel bit,
432 int verb) const {
433 MoFEM::Interface &m_field = cOre;
434 auto ref_ents_ptr = m_field.get_ref_ents();
435 auto ref_fe_ptr = m_field.get_ref_finite_elements();
437
438 for (Range::const_pair_iterator pit = ents.pair_begin();
439 pit != ents.pair_end(); pit++) {
440 // get first and last element of range
441 EntityHandle f = pit->first;
442 EntityHandle s = pit->second;
443 Range seed_ents_range; // entities seeded not in database
444 // find ents to add
445 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
446 .findEntsToAdd(f, s, seed_ents_range);
447 // add elements
448 if (!seed_ents_range.empty())
449 CHKERR SetBitRefLevelTool(m_field, bit, ref_ents_ptr, ref_fe_ptr)
450 .addEntsToDatabase(seed_ents_range);
451 }
453}
454
456 const std::string field_name, const BitRefLevel bit, int verb) const {
457 MoFEM::Interface &m_field = cOre;
459 EntityHandle field_meshset = m_field.get_field_meshset(field_name);
460 Range field_ents;
461 CHKERR m_field.get_moab().get_entities_by_handle(field_meshset, field_ents,
462 true);
463 CHKERR setEntitiesBitRefLevel(field_ents, bit, verb);
465}
466
468 const EntityType type, const BitRefLevel bit, const BitRefLevel mask,
469 int verb) const {
471 Range ents;
473 CHKERR setBitRefLevel(ents, BitRefLevel(), false, verb);
475}
476
478 const int dim, const BitRefLevel bit, const BitRefLevel mask,
479 int verb) const {
481 Range ents;
483 CHKERR setBitRefLevel(ents, BitRefLevel(), false, verb);
485}
486
488 const BitRefLevel bit,
489 int verb) const {
490 MoFEM::Interface &m_field = cOre;
491 auto ref_ents_ptr = m_field.get_ref_ents();
492 auto ref_fe_ptr = m_field.get_ref_finite_elements();
494 // Add ref entity
495 std::pair<RefEntity_multiIndex::iterator, bool> p_ent =
496 const_cast<RefEntity_multiIndex *>(ref_ents_ptr)
497 ->insert(boost::shared_ptr<RefEntity>(
498 new RefEntity(m_field.get_basic_entity_data_ptr(), meshset)));
499 *(const_cast<RefEntity *>(p_ent.first->get())->getBitRefLevelPtr()) |= bit;
500 // Add ref element
501 boost::shared_ptr<RefElement> fe_ptr =
502 boost::shared_ptr<RefElement>(new RefElement_MESHSET(*p_ent.first));
503 std::pair<RefElement_multiIndex::iterator, bool> p_fe =
504 const_cast<RefElement_multiIndex *>(ref_fe_ptr)->insert(fe_ptr);
505
507 MOFEM_LOG("BitRefSelf", Sev::noisy)
508 << "Add meshset as ref_ent " << **p_fe.first;
509
511}
512
514 const int dim,
515 const BitRefLevel bit,
516 int verb) const {
517 MoFEM::Interface &m_field = cOre;
519 Range ents;
520 CHKERR m_field.get_moab().get_entities_by_dimension(meshset, dim, ents,
521 false);
522 CHKERR setBitRefLevel(ents, bit, false, verb);
524}
525
527 const EntityType type,
528 const BitRefLevel bit,
529 int verb) const {
530 MoFEM::Interface &m_field = cOre;
532 Range ents;
533 CHKERR m_field.get_moab().get_entities_by_type(meshset, type, ents, false);
534 CHKERR setBitRefLevel(ents, bit, false, verb);
536}
537
539 const BitRefLevel bit,
540 int verb) const {
541 MoFEM::Interface &m_field = cOre;
543 std::vector<const BitRefLevel *> ents_bits_vec;
544 CHKERR RefEntity::getBitRefLevel(m_field.get_moab(), ents, ents_bits_vec);
545 for (auto it : ents_bits_vec)
546 const_cast<BitRefLevel &>(*it) |= bit;
548}
549
551 const int dim,
552 const BitRefLevel bit,
553 int verb) const {
554 MoFEM::Interface &m_field = cOre;
555 moab::Interface &moab = m_field.get_moab();
556 Range ents, adj;
558 CHKERR moab.get_entities_by_dimension(meshset, dim, ents, true);
559 for (int dd = dim - 1; dd >= 0; dd--) {
560 CHKERR moab.get_adjacencies(ents, dd, false, adj, moab::Interface::UNION);
561 }
562 ents.merge(adj);
563 if (verb == VERY_NOISY) {
565 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Add add bit ref level by dim ";
566 }
567
568 CHKERR addBitRefLevel(ents, bit, verb);
570}
571
573 const bool b, int verb) const {
574 MoFEM::Interface &m_field = cOre;
576 std::vector<const BitRefLevel *> ents_bits_vec;
577 CHKERR RefEntity::getBitRefLevel(m_field.get_moab(), ents, ents_bits_vec);
578 for (auto it = ents_bits_vec.begin(); it != ents_bits_vec.end(); ++it) {
579 const_cast<BitRefLevel &>(**it)[n] = b;
580 }
581 if (verb == VERY_NOISY) {
582 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Set bit to " << ents;
583 }
585}
586
588 int verb) const {
589 MoFEM::Interface &m_field = cOre;
590 auto ref_ents_ptr = m_field.get_ref_ents();
592 auto hi_dit = ref_ents_ptr->end();
593 for (auto dit = ref_ents_ptr->begin(); dit != hi_dit; ++dit) {
594 (*const_cast<RefEntity *>(dit->get())->getBitRefLevelPtr())[n] = b;
595 if (verb >= VERY_VERBOSE) {
596 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Set bit to " << **dit;
597 }
598 }
600}
601
603 const BitRefLevel mask,
604 int verb) const {
606 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
608}
609
611 const BitRefLevel mask, int verb,
612 MoFEMTypes mf) const {
613 MoFEM::Interface &m_field = cOre;
614 auto ref_ents_ptr = m_field.get_ref_ents();
616 RefEntity_change_right_shift right_shift(1, mask);
617 for (int ii = 0; ii < shift; ii++) {
618 // delete bits on the right which are shifted to zero
619 BitRefLevel delete_bits = BitRefLevel().set(0) & mask;
620 if (delete_bits.any()) {
621 CHKERR m_field.delete_ents_by_bit_ref(delete_bits, delete_bits, true,
622 verb, mf);
623 }
624 for (RefEntity_multiIndex::iterator ent_it = ref_ents_ptr->begin();
625 ent_it != ref_ents_ptr->end(); ent_it++) {
626 if (verb >= NOISY) {
628 MOFEM_LOG("BitRefSelf", Sev::noisy)
629 << (*ent_it)->getBitRefLevel() << " : ";
630 }
631 right_shift(const_cast<boost::shared_ptr<RefEntity> &>(*ent_it));
632 if (verb == VERY_NOISY) {
634 MOFEM_LOG("BitRefSelf", Sev::noisy) << (*ent_it)->getBitRefLevel();
635 }
636 }
637 }
639}
640
642 const BitRefLevel mask,
643 const char *file_name,
644 const char *file_type,
645 const char *options,
646 const bool check_for_empty) const {
647 MoFEM::Interface &m_field = cOre;
648 moab::Interface &moab(m_field.get_moab());
650 EntityHandle meshset;
651 CHKERR moab.create_meshset(MESHSET_SET, meshset);
652 CHKERR getEntitiesByRefLevel(bit, mask, meshset);
653 int nb_ents;
654 CHKERR moab.get_number_entities_by_handle(meshset, nb_ents, true);
655 if (check_for_empty && !nb_ents) {
656 MOFEM_LOG("SELF", Sev::warning)
657 << "No entities to save < " << file_name << " > in writeBitLevel";
659 }
660
661 CHKERR moab.write_file(file_name, file_type, options, &meshset, 1);
662 CHKERR moab.delete_entities(&meshset, 1);
664}
665
668 const int dim, const char *file_name,
669 const char *file_type, const char *options,
670 const bool check_for_empty) const {
671 MoFEM::Interface &m_field = cOre;
672 moab::Interface &moab(m_field.get_moab());
674 Range ents;
676 if (check_for_empty && ents.empty()) {
677 MOFEM_LOG("SELF", Sev::warning)
678 << "No entities to save < " << file_name << " > in writeBitLevelByDim";
680 }
681 EntityHandle meshset;
682 CHKERR moab.create_meshset(MESHSET_SET, meshset);
683 CHKERR moab.add_entities(meshset, ents);
684 CHKERR moab.write_file(file_name, file_type, options, &meshset, 1);
685 CHKERR moab.delete_entities(&meshset, 1);
687}
688
690 const BitRefLevel bit, const BitRefLevel mask, const EntityType type,
691 const char *file_name, const char *file_type, const char *options,
692 const bool check_for_empty) const {
693 MoFEM::Interface &m_field = cOre;
694 moab::Interface &moab(m_field.get_moab());
696 Range ents;
698 if (check_for_empty && ents.empty()) {
699 MOFEM_LOG("SELF", Sev::warning)
700 << "No entities to save < " << file_name << " > in writeBitLevelByType";
702 }
703 EntityHandle meshset;
704 CHKERR moab.create_meshset(MESHSET_SET, meshset);
705 CHKERR moab.add_entities(meshset, ents);
706 CHKERR moab.write_file(file_name, file_type, options, &meshset, 1);
707 CHKERR moab.delete_entities(&meshset, 1);
709}
710
712 const char *file_name, const char *file_type, const char *options,
713 const bool check_for_empty) const {
714 MoFEM::Interface &m_field = cOre;
715 moab::Interface &moab(m_field.get_moab());
717 EntityHandle meshset;
718 Range ents;
720 if (check_for_empty && ents.empty())
722 CHKERR moab.create_meshset(MESHSET_SET, meshset);
723 CHKERR moab.add_entities(meshset, ents);
724 CHKERR moab.write_file(file_name, file_type, options, &meshset, 1);
725 CHKERR moab.delete_entities(&meshset, 1);
727}
728
730 const BitRefLevel mask, const EntityType type, const char *file_name,
731 const char *file_type, const char *options) {
733 for (int ll = 0; ll != BITREFLEVEL_SIZE; ++ll) {
734 std::string name = boost::lexical_cast<std::string>(ll) + "_" + file_name;
735 CHKERR writeBitLevelByType(BitRefLevel().set(ll), mask, type, name.c_str(),
736 file_type, options, true);
737 }
739}
740
742 const BitRefLevel bit, const BitRefLevel mask, const EntityType type,
743 const EntityHandle meshset, int verb) const {
744 MoFEM::Interface &m_field = cOre;
745 moab::Interface &moab(m_field.get_moab());
747 Range ents;
748 CHKERR getEntitiesByTypeAndRefLevel(bit, mask, type, ents, verb);
749 CHKERR moab.add_entities(meshset, ents);
751}
752
754 const BitRefLevel mask,
755 Range &ents,
756 int verb) const {
757 MoFEM::Interface &m_field = cOre;
758 moab::Interface &moab(m_field.get_moab());
760
761 std::vector<EntityHandle> ents_vec;
762 ents_vec.reserve(ents.size());
763
764 std::vector<BitRefLevel *> tags_bits_ptr_vec(ents.size());
765
766 Range swap_ents;
767 auto hint = swap_ents.begin();
768
769 for (Range::pair_iterator p_eit = ents.pair_begin(); p_eit != ents.pair_end();
770 ++p_eit) {
771
772 EntityHandle f = p_eit->first;
773 const EntityHandle s = p_eit->second;
774
775 // get bits on entities
776 rval = moab.tag_get_by_ptr(cOre.get_th_RefBitLevel(), Range(f, s),
777 (const void **)(&*tags_bits_ptr_vec.begin()));
778
779 if (rval == MB_SUCCESS) {
780
781 auto bit_it = tags_bits_ptr_vec.begin();
782
783 auto check = [&bit, &mask](const auto &entity_bit) -> bool {
784 return
785
786 (entity_bit & bit).any() &&
787
788 ((entity_bit & mask) == entity_bit);
789 };
790
791 while (f != s + 1) {
792
793 while (f != s + 1 && !check(**bit_it)) {
794 ++bit_it;
795 ++f;
796 }
797
798 if (f != s + 1) {
799
800 const EntityHandle start = f;
801
802 while (f != (s + 1) && check(**bit_it)) {
803 ++bit_it;
804 ++f;
805 };
806
807 hint = swap_ents.insert(hint, start, f - 1);
808 }
809 }
810 }
811 }
812
813 ents.swap(swap_ents);
814
816}
817
819 const BitRefLevel bit, const BitRefLevel mask, const EntityType type,
820 Range &ents, int verb) const {
821 MoFEM::Interface &m_field = cOre;
822 moab::Interface &moab(m_field.get_moab());
824 CHKERR moab.get_entities_by_type(0, type, ents, false);
825 CHKERR filterEntitiesByRefLevel(bit, mask, ents, verb);
827}
828
830 const BitRefLevel bit, const BitRefLevel mask, const int dim,
831 const EntityHandle meshset, int verb) const {
832 MoFEM::Interface &m_field = cOre;
833 moab::Interface &moab(m_field.get_moab());
835 Range ents;
836 CHKERR getEntitiesByDimAndRefLevel(bit, mask, dim, ents, verb);
837 CHKERR moab.add_entities(meshset, ents);
839}
840
842 const BitRefLevel bit, const BitRefLevel mask, const int dim, Range &ents,
843 int verb) const {
844 MoFEM::Interface &m_field = cOre;
845 moab::Interface &moab(m_field.get_moab());
847 CHKERR moab.get_entities_by_dimension(0, dim, ents, false);
848 CHKERR filterEntitiesByRefLevel(bit, mask, ents, verb);
850}
851
853 const BitRefLevel mask,
854 const EntityHandle meshset,
855 const int verb) const {
856 MoFEM::Interface &m_field = cOre;
857 moab::Interface &moab(m_field.get_moab());
859 Range ents;
860 CHKERR getEntitiesByRefLevel(bit, mask, ents, verb);
861 CHKERR moab.add_entities(meshset, ents);
863}
864
866 const BitRefLevel mask,
867 Range &ents,
868 const int verb) const {
869 MoFEM::Interface &m_field = cOre;
870 moab::Interface &moab(m_field.get_moab());
872 Range meshset_ents;
873 CHKERR moab.get_entities_by_type(0, MBENTITYSET, meshset_ents, false);
874 CHKERR moab.get_entities_by_handle(0, ents, false);
875 ents.merge(meshset_ents);
876 CHKERR filterEntitiesByRefLevel(bit, mask, ents, verb);
878}
879
881 const BitRefLevel mask,
882 const EntityType type,
883 Range &ents,
884 const int verb) const {
885 MoFEM::Interface &m_field = cOre;
886 auto ref_ents_ptr = m_field.get_ref_ents();
888 auto &ref_ents = ref_ents_ptr->get<Ent_Ent_mi_tag>();
889 auto it = ref_ents.lower_bound(get_id_for_min_type(type));
890 auto hi_it = ref_ents.upper_bound(get_id_for_max_type(type));
891 std::vector<EntityHandle> ents_vec;
892 ents_vec.reserve(std::distance(it, hi_it));
893 for (; it != hi_it; it++) {
894 const BitRefLevel &ent_bit = it->get()->getBitRefLevel();
895 if ((ent_bit & mask) == ent_bit && (ent_bit & bit).any())
896 ents_vec.emplace_back(it->get()->getEnt());
897 }
898 ents.insert_list(ents_vec.begin(), ents_vec.end());
899 if (verb > NOISY)
900 MOFEM_LOG("BitRefSelf", Sev::noisy)
901 << "getEntitiesByParentType: " << ents << endl;
903}
904
906 MoFEM::Interface &m_field = cOre;
907 moab::Interface &moab = m_field.get_moab();
909 CHKERR moab.get_entities_by_handle(0, ents, false);
910 ents = subtract(ents, ents.subset_by_type(MBENTITYSET));
913}
914
916 MoFEM::Interface &m_field = cOre;
917 auto ref_ents_ptr = m_field.get_ref_ents();
919 auto eit = ents.begin();
920 for (; eit != ents.end();) {
921 auto rit = ref_ents_ptr->get<Ent_mi_tag>().find(*eit);
922 if (rit != ref_ents_ptr->get<Ent_mi_tag>().end()) {
923 eit = ents.erase(eit);
924 } else {
925 eit++;
926 }
927 }
929}
930
932BitRefManager::getAdjacenciesEquality(const EntityHandle from_entity,
933 const int to_dimension,
934 Range &adj_entities) const {
935 MoFEM::Interface &m_field = cOre;
936 moab::Interface &moab(m_field.get_moab());
938 BitRefLevel bit_from_entity;
939 CHKERR moab.tag_get_data(cOre.get_th_RefBitLevel(), &from_entity, 1,
940 &bit_from_entity);
941 CHKERR moab.get_adjacencies(&from_entity, 1, to_dimension, false,
942 adj_entities);
943 std::vector<BitRefLevel> bit_levels(adj_entities.size());
944 CHKERR moab.tag_get_data(cOre.get_th_RefBitLevel(), adj_entities,
945 &*bit_levels.begin());
946 auto b_it = bit_levels.begin();
947 auto eit = adj_entities.begin();
948 for (; eit != adj_entities.end(); b_it++) {
949 if (bit_from_entity != *b_it) {
950 eit = adj_entities.erase(eit);
951 } else {
952 eit++;
953 }
954 }
956}
957
959 const int to_dimension,
960 Range &adj_entities) const {
961 MoFEM::Interface &m_field = cOre;
962 moab::Interface &moab(m_field.get_moab());
964 BitRefLevel bit_from_entity;
965 CHKERR moab.tag_get_data(cOre.get_th_RefBitLevel(), &from_entity, 1,
966 &bit_from_entity);
967 CHKERR moab.get_adjacencies(&from_entity, 1, to_dimension, false,
968 adj_entities);
969 std::vector<BitRefLevel> bit_levels(adj_entities.size());
970 CHKERR moab.tag_get_data(cOre.get_th_RefBitLevel(), adj_entities,
971 &*bit_levels.begin());
972 std::vector<BitRefLevel>::iterator b_it = bit_levels.begin();
973 Range::iterator eit = adj_entities.begin();
974 for (; eit != adj_entities.end(); b_it++) {
975 if (!(bit_from_entity & (*b_it)).any()) {
976 eit = adj_entities.erase(eit);
977 } else {
978 eit++;
979 }
980 }
982}
983
985 const Problem *problem_ptr, const EntityHandle *from_entities,
986 const int num_entities, const int to_dimension, Range &adj_entities,
987 const int operation_type, const int verb) const {
989 BitRefLevel bit = problem_ptr->getBitRefLevel();
990 CHKERR getAdjacencies(bit, from_entities, num_entities, to_dimension,
991 adj_entities, operation_type);
993}
994
996 const BitRefLevel bit, const EntityHandle *from_entities,
997 const int num_entities, const int to_dimension, Range &adj_entities,
998 const int operation_type, const int verb) const {
999 MoFEM::Interface &m_field = cOre;
1000 moab::Interface &moab(m_field.get_moab());
1002
1003 if (verb > VERBOSE) {
1005 MOFEM_LOG("BitRef", Sev::noisy) << "from: " << bit;
1006 }
1007
1008 CHKERR moab.get_adjacencies(from_entities, num_entities, to_dimension, false,
1009 adj_entities, operation_type);
1010 std::vector<BitRefLevel> bit_levels(adj_entities.size());
1011 CHKERR moab.tag_get_data(cOre.get_th_RefBitLevel(), adj_entities,
1012 &*bit_levels.begin());
1013 std::vector<BitRefLevel>::iterator b_it = bit_levels.begin();
1014 for (Range::iterator eit = adj_entities.begin(); eit != adj_entities.end();
1015 b_it++) {
1016
1017 if (verb > VERBOSE) {
1018 RefEntity adj_entity(m_field.get_basic_entity_data_ptr(), *eit);
1020 MOFEM_LOG("BitRef", Sev::noisy)
1021 << "to: " << adj_entity.getBitRefLevel() << " : " << adj_entity;
1022 }
1023
1024 if (!((*b_it) & bit).any()) {
1025 eit = adj_entities.erase(eit);
1026 } else {
1027 eit++;
1028 }
1029 }
1030 if (b_it != bit_levels.end()) {
1031 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Data inconsistency");
1032 }
1034}
1035
1037 const EntityHandle parent, const BitRefLevel &parent_bit,
1038 const BitRefLevel &parent_mask, const BitRefLevel &child_bit,
1039 const BitRefLevel &child_mask, const EntityHandle child,
1040 EntityType child_type, const bool recursive, int verb) {
1041 MoFEM::Interface &m_field = cOre;
1042 moab::Interface &moab = m_field.get_moab();
1044 Range parent_ents;
1045 CHKERR moab.get_entities_by_handle(parent, parent_ents, recursive);
1046 CHKERR filterEntitiesByRefLevel(parent_bit, parent_mask, parent_ents, verb);
1047 if (verb >= VERY_VERBOSE) {
1049 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Parnets: " << parent;
1050 }
1051 Range children_ents;
1052 CHKERR updateRangeByChildren(parent_ents, children_ents);
1053 if (child_type < MBMAXTYPE)
1054 children_ents = children_ents.subset_by_type(child_type);
1055 CHKERR filterEntitiesByRefLevel(child_bit, BitRefLevel().set(), children_ents,
1056 verb);
1057 CHKERR moab.add_entities(child, children_ents);
1059}
1060
1062 const EntityHandle parent, const BitRefLevel &child_bit,
1063 const EntityHandle child, EntityType child_type, const bool recursive,
1064 int verb) {
1067 parent, BitRefLevel().set(), BitRefLevel().set(), child_bit, child_bit,
1068 child, child_type, recursive, verb);
1070}
1071
1073 const BitRefLevel &child_bit, int verb) {
1074 MoFEM::Interface &m_field = cOre;
1075 moab::Interface &moab = m_field.get_moab();
1076 auto fields_ptr = m_field.get_fields();
1077 auto ref_ents_ptr = m_field.get_ref_ents();
1079
1080 for (auto &fit : (*fields_ptr)) {
1081
1082 const EntityHandle meshset = fit->getMeshset();
1083 Range parent_ents;
1084 CHKERR moab.get_entities_by_handle(meshset, parent_ents, true);
1085
1086 if (verb >= VERY_VERBOSE) {
1088 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Parnets:" << std::endl
1089 << parent_ents << std::endl;
1090 }
1091
1092 Range children_ents;
1093 CHKERR updateRangeByChildren(parent_ents, children_ents);
1094 CHKERR filterEntitiesByRefLevel(child_bit, BitRefLevel().set(),
1095 children_ents, verb);
1096
1097 CHKERR moab.add_entities(meshset, children_ents);
1098
1099 if (verb >= VERY_VERBOSE) {
1101 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Children: " << std::endl
1102 << children_ents << std::endl;
1103 }
1104 }
1106}
1107
1109 const std::string name, const BitRefLevel &child_bit, int verb) {
1110 MoFEM::Interface &m_field = cOre;
1111 moab::Interface &moab = m_field.get_moab();
1113
1114 EntityHandle field_meshset = m_field.get_field_structure(name)->getMeshset();
1115
1116 Range parent_ents;
1117 CHKERR moab.get_entities_by_handle(field_meshset, parent_ents, true);
1118
1119 if (verb >= VERBOSE) {
1121 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Parnets:" << endl
1122 << parent_ents << std::endl;
1123 }
1124
1125 Range children_ents;
1126 CHKERR updateRangeByChildren(parent_ents, children_ents);
1127 CHKERR filterEntitiesByRefLevel(child_bit, BitRefLevel().set(), children_ents,
1128 verb);
1129
1130 CHKERR moab.add_entities(field_meshset, children_ents);
1131
1132 if (verb >= VERBOSE) {
1134 MOFEM_LOG("BitRefSelf", Sev::noisy) << "Children: " << endl
1135 << children_ents << std::endl;
1136 }
1137
1139}
1140
1142 const std::string name, const BitRefLevel &child_bit,
1143 const EntityType fe_ent_type, int verb) {
1144 MoFEM::Interface &m_field = cOre;
1146 EntityHandle meshset = m_field.get_finite_element_meshset(name);
1147 CHKERR updateMeshsetByEntitiesChildren(meshset, child_bit, meshset,
1148 fe_ent_type, false, verb);
1150}
1151
1153 Range &child_ents, MoFEMTypes bh) {
1154 MoFEM::Interface &m_field = cOre;
1155 auto ref_ents_ptr = m_field.get_ref_ents();
1157 auto &ref_ents =
1158 const_cast<RefEntity_multiIndex *>(ref_ents_ptr)->get<Ent_Ent_mi_tag>();
1159 std::vector<EntityHandle> child_ents_vec;
1160 child_ents_vec.reserve(ref_ents.size());
1161 for (Range::const_pair_iterator pit = parent_ents.const_pair_begin();
1162 pit != parent_ents.const_pair_end(); pit++) {
1163 const auto f = pit->first;
1164 auto it = ref_ents.lower_bound(f);
1165 if (it != ref_ents.end()) {
1166 const auto s = pit->second;
1167 auto hi_it = ref_ents.upper_bound(s);
1168 if (bh == MF_EXIST) {
1169 if (std::distance(it, hi_it) != (s - f) + 1) {
1170 SETERRQ2(
1171 PETSC_COMM_SELF, MOFEM_NOT_FOUND,
1172 "Number of entities and entities parents is different %d != %d ",
1173 std::distance(it, hi_it), (s - f) + 1);
1174 }
1175 }
1176 for (; it != hi_it; it++) {
1177#ifndef NDEBUG
1178 if (it->get()->getEntType() == MBENTITYSET) {
1179 SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSIBLE_CASE,
1180 "This should not happen; Entity should not have part of the "
1181 "meshset. It has no children.");
1182 }
1183#endif
1184 child_ents_vec.emplace_back((*it)->getEnt());
1185 }
1186 } else if (bh == MF_EXIST) {
1187 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_FOUND, "Entities not found");
1188 }
1189 }
1190 child_ents.insert_list(child_ents_vec.begin(), child_ents_vec.end());
1192}
1193
1195 Range &parent_ents,
1196 MoFEMTypes bh) {
1197 MoFEM::Interface &m_field = cOre;
1198 auto ref_ents_ptr = m_field.get_ref_ents();
1200 auto &ref_ents =
1201 const_cast<RefEntity_multiIndex *>(ref_ents_ptr)->get<Ent_mi_tag>();
1202 std::vector<EntityHandle> parent_ents_vec;
1203 parent_ents_vec.reserve(ref_ents.size());
1204 for (Range::const_pair_iterator pit = child_ents.const_pair_begin();
1205 pit != child_ents.const_pair_end(); pit++) {
1206 const auto f = pit->first;
1207 auto it = ref_ents.lower_bound(f);
1208 if (it != ref_ents.end()) {
1209 const auto s = pit->second;
1210 auto hi_it = ref_ents.upper_bound(s);
1211 if (bh == MF_EXIST) {
1212 if (std::distance(it, hi_it) != (s - f) + 1) {
1213 SETERRQ2(
1214 PETSC_COMM_SELF, MOFEM_NOT_FOUND,
1215 "Number of entities and enties parents is diffrent %d != %d ",
1216 std::distance(it, hi_it), (s - f) + 1);
1217 }
1218 }
1219 for (; it != hi_it; it++) {
1220#ifndef NDEBUG
1221 if (it->get()->getEntType() == MBENTITYSET) {
1222 SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSIBLE_CASE,
1223 "This should not happen; Entity should not have part of the "
1224 "meshset. It has no children.");
1225 }
1226#endif
1227 auto parent = (*it)->getParentEnt();
1228 if (parent)
1229 parent_ents_vec.emplace_back(parent);
1230 }
1231 } else if (bh == MF_EXIST) {
1232 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_FOUND, "Entities not found");
1233 }
1234 }
1235 parent_ents.insert_list(parent_ents_vec.begin(), parent_ents_vec.end());
1237}
1238
1239} // namespace MoFEM
static Index< 'n', 3 > n
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:314
@ VERY_VERBOSE
Definition: definitions.h:223
@ QUIET
Definition: definitions.h:221
@ NOISY
Definition: definitions.h:224
@ VERBOSE
Definition: definitions.h:222
@ VERY_NOISY
Definition: definitions.h:225
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
Definition: definitions.h:110
@ MF_EXIST
Definition: definitions.h:113
#define BITREFLEVEL_SIZE
max number of refinements
Definition: definitions.h:232
#define MOAB_THROW(err)
Check error code of MoAB function and throw MoFEM exception.
Definition: definitions.h:554
#define MAX_CORE_TMP
maximal number of cores
Definition: definitions.h:230
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
@ MOFEM_NOT_FOUND
Definition: definitions.h:46
@ MOFEM_IMPOSIBLE_CASE
Definition: definitions.h:48
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:44
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:45
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
const int dim
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
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
virtual const RefEntity_multiIndex * get_ref_ents() const =0
Get the ref ents object.
virtual const Field_multiIndex * get_fields() const =0
Get the fields object.
virtual const RefElement_multiIndex * get_ref_finite_elements() const =0
Get the ref finite elements object.
virtual MoFEMErrorCode getAdjacenciesAny(const EntityHandle from_entity, const int to_dimension, Range &adj_entities) const
Get the adjacencies associated with a entity to entities of a specified dimension.
MoFEMErrorCode updateMeshsetByEntitiesChildren(const EntityHandle parent, const BitRefLevel &parent_bit, const BitRefLevel &parent_mask, const BitRefLevel &child_bit, const BitRefLevel &child_mask, const EntityHandle child, EntityType child_type, const bool recursive=false, int verb=0)
Get child entities form meshset containing parent entities.
MoFEMErrorCode shiftLeftBitRef(const int shift, const BitRefLevel mask=BitRefLevel().set(), int verb=DEFAULT_VERBOSITY) const
left shift bit ref levelthis results of deletion of entities on far left side
MoFEMErrorCode setElementsBitRefLevel(const Range &ents, const BitRefLevel bit=BitRefLevel(), int verb=QUIET) const
add entities to database and set bit ref level
MoFEMErrorCode shiftRightBitRef(const int shift, const BitRefLevel mask=BitRefLevel().set(), int verb=DEFAULT_VERBOSITY, MoFEMTypes mf=MF_ZERO) const
right shift bit ref level
MoFEMErrorCode setEntitiesBitRefLevel(const Range &ents, const BitRefLevel bit=BitRefLevel(), int verb=QUIET) const
add entities to database and set bit ref level
MoFEMErrorCode getEntitiesByDimAndRefLevel(const BitRefLevel bit, const BitRefLevel mask, const int dim, const EntityHandle meshset, int verb=0) const
add all ents from ref level given by bit to meshset
MoFEMErrorCode setBitRefLevel(const Range &ents, const BitRefLevel bit, const bool only_tets=true, int verb=0, Range *adj_ents_ptr=nullptr) const
add entities to database and set bit ref level
MoFEMErrorCode setBitRefLevelByType(const EntityHandle meshset, const EntityType type, const BitRefLevel bit, int verb=QUIET) const
Set the Bit Ref Level By Type object.
MoFEMErrorCode addBitRefLevelByDim(const EntityHandle meshset, const int dim, const BitRefLevel bit, int verb=QUIET) const
add bit ref level by dimension
MoFEMErrorCode setBitLevelToMeshset(const EntityHandle meshset, const BitRefLevel bit, int verb=0) const
MoFEMErrorCode setNthBitRefLevel(const Range &ents, const int n, const bool b, int verb=QUIET) const
Set nth bit ref level.
MoFEMErrorCode addToDatabaseBitRefLevelByType(const EntityType type, const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set(), int verb=QUIET) const
Add entities which exist in MoAB database, and have set appropiate BitRef level tag,...
MoFEMErrorCode updateFiniteElementMeshsetByEntitiesChildren(const std::string name, const BitRefLevel &child_bit, const EntityType fe_ent_type, int verb=0)
update finite element meshset by child entities
MoFEMErrorCode addToDatabaseBitRefLevelByDim(const int dim, const BitRefLevel bit, const BitRefLevel mask=BitRefLevel().set(), int verb=QUIET) const
Add entities which exist in MoAB database, and have set appropiate BitRef level tag,...
virtual MoFEMErrorCode getAdjacencies(const Problem *problem_ptr, const EntityHandle *from_entities, const int num_entities, const int to_dimension, Range &adj_entities, const int operation_type=moab::Interface::INTERSECT, const int verb=0) const
Get the adjacencies associated with a entity to entities of a specified dimension.
MoFEMErrorCode setFieldEntitiesBitRefLevel(const std::string field_name, const BitRefLevel bit=BitRefLevel(), int verb=QUIET) const
Set the bit ref level to entities in the field meshset.
MoFEMErrorCode getEntitiesByTypeAndRefLevel(const BitRefLevel bit, const BitRefLevel mask, const EntityType type, const EntityHandle meshset, int verb=0) const
add all ents from ref level given by bit to meshset
MoFEMErrorCode addBitRefLevel(const Range &ents, const BitRefLevel bit, int verb=QUIET) const
add bit ref level to ref entity
MoFEMErrorCode getEntitiesByRefLevel(const BitRefLevel bit, const BitRefLevel mask, const EntityHandle meshset, const int verb=QUIET) const
add all ents from ref level given by bit to meshset
MoFEMErrorCode updateFieldMeshsetByEntitiesChildren(const BitRefLevel &child_bit, int verb=0)
update fields meshesets by child entities
virtual MoFEMErrorCode getAdjacenciesEquality(const EntityHandle from_entity, const int to_dimension, Range &adj_entities) const
Get the adjacencies associated with a entity to entities of a specified dimension.
MoFEMErrorCode updateRangeByParent(const Range &parent, Range &child, MoFEMTypes bh=MF_ZERO)
Update range by prents.
MoFEMErrorCode setBitRefLevelByDim(const EntityHandle meshset, const int dim, const BitRefLevel bit, int verb=QUIET) const
Set the Bit Ref Level By Dim object.
MoFEMErrorCode updateRangeByChildren(const Range &parent, Range &child, MoFEMTypes bh=MF_ZERO)
Update range by childresn.
virtual EntityHandle get_finite_element_meshset(const std::string &name) const =0
virtual Field * get_field_structure(const std::string &name)=0
get field structure
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
Definition: LogManager.cpp:364
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:311
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:342
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
Definition: LogManager.hpp:287
#define MOFEM_LOG_ATTRIBUTES(channel, bit)
Add attributes to channel.
Definition: LogManager.hpp:299
#define MOFEM_LOG_FUNCTION()
Set scope.
Definition: LogManager.hpp:328
auto bit
set bit
const Tensor1_Expr< const dTensor0< T, Dim, i >, typename promote< T, double >::V, Dim, i > d(const Tensor0< T * > &a, const Index< i, Dim > index, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: dTensor0.hpp:27
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Exceptions.hpp:85
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
EntityHandle get_id_for_max_type()
EntityHandle get_id_for_min_type()
moab::Range::iterator insertOrdered(Range &r, Extractor, Iterator begin_iter, Iterator end_iter)
Insert ordered mofem multi-index into range.
Definition: Templates.hpp:1410
RefEntityTmp< 0 > RefEntity
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1965
constexpr auto field_name
Managing BitRefLevels.
MoFEMErrorCode getAllEntitiesNotInDatabase(Range &ents) const
Get all entities not in database.
MoFEMErrorCode writeBitLevelByDim(const BitRefLevel bit, const BitRefLevel mask, const int dim, const char *file_name, const char *file_type, const char *options, const bool check_for_empty=true) const
Write bit ref level to file.
MoFEMErrorCode writeBitLevelByType(const BitRefLevel bit, const BitRefLevel mask, const EntityType type, const char *file_name, const char *file_type, const char *options, const bool check_for_empty=true) const
Write bit ref level to file.
MoFEMErrorCode getEntitiesByParentType(const BitRefLevel bit, const BitRefLevel mask, const EntityType type, Range &ents, const int verb=QUIET) const
get entities by bit ref level and type of parent
MoFEMErrorCode filterEntitiesByRefLevel(const BitRefLevel bit, const BitRefLevel mask, Range &ents, int verb=QUIET) const
filter entities by bit ref level
BitRefManager(const MoFEM::Core &core)
MoFEMErrorCode filterEntitiesNotInDatabase(Range &ents) const
Get entities not in database.
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
MoFEMErrorCode writeEntitiesAllBitLevelsByType(const BitRefLevel mask, const EntityType type, const char *file_name, const char *file_type, const char *options)
Write all entities by bit levels and type.
MoFEMErrorCode writeBitLevel(const BitRefLevel bit, const BitRefLevel mask, const char *file_name, const char *file_type, const char *options, const bool check_for_empty=true) const
Write bit ref level to file.
MoFEMErrorCode writeEntitiesNotInDatabase(const char *file_name, const char *file_type, const char *options, const bool check_for_empty=true) const
Write ents not in database.
virtual moab::Interface & get_moab()=0
virtual EntityHandle get_field_meshset(const std::string name) const =0
get field meshset
virtual const int getValue() const =0
Get the core.
virtual MoFEMErrorCode delete_ents_by_bit_ref(const BitRefLevel bit, const BitRefLevel mask, const bool remove_parent=false, int verb=DEFAULT_VERBOSITY, MoFEMTypes mf=MF_ZERO)=0
delete entities form mofem and moab database
virtual boost::shared_ptr< BasicEntityData > & get_basic_entity_data_ptr()=0
Get pointer to basic entity data.
Core (interface) class.
Definition: Core.hpp:92
Tag get_th_RefBitLevel() const
Definition: Core.hpp:208
Deprecated interface functions.
EntityHandle getMeshset() const
Get field meshset.
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
Definition: LogManager.cpp:279
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
Definition: LogManager.cpp:323
static boost::shared_ptr< std::ostream > getStrmSync()
Get the strm sync object.
Definition: LogManager.cpp:327
static bool checkIfChannelExist(const std::string channel)
Check if channel exist.
Definition: LogManager.cpp:374
static boost::shared_ptr< std::ostream > getStrmSelf()
Get the strm self object.
Definition: LogManager.cpp:319
keeps basic data about problem
BitRefLevel getBitRefLevel() const
keeps data about abstract EDGE finite element
keeps data about abstract MESHSET finite element
keeps data about abstract PRISM finite element
keeps data about abstract VERTEX finite element
keeps data about abstract TRI finite element
keeps data about abstract TET finite element
Extract entity handle form multi-index container.
Definition: Templates.hpp:1386
Struct keeps handle to refined handle.
const BitRefLevel & getBitRefLevel() const
Get entity ref bit refinement signature.
static MoFEMErrorCode getBitRefLevel(Interface &moab, Range ents, std::vector< BitRefLevel > &vec_bit_ref_level)
tool class with methods used more than twp times
boost::shared_ptr< BasicEntityData > & baseEntData
base entity data
MoFEM::Interface & mField
SetBitRefLevelTool(MoFEM::Interface &m_field, const BitRefLevel &bit, const RefEntity_multiIndex *ref_ents_ptr, const RefElement_multiIndex *ref_element_ptr)
constrictor
const RefEntity_multiIndex * refEntsPtr
access to ents database
MoFEMErrorCode findEntsToAdd(EntityHandle f, EntityHandle s, Range &seed_ents_range) const
find entities and change entity bit if in database
MoFEMErrorCode addEntsToDatabase(const Range &seed_ents_range) const
const BitRefLevel & bIt
bit to set
MoFEMErrorCode addElementsToDatabase(Range &seed_fe_range) const
const RefElement_multiIndex * refElementPtr
access to fe database
MoFEMErrorCode addEntsToDatabaseImpl(const Range &seed_ents_range) const
add entities to database
MoFEMErrorCode findElementsToAdd(EntityHandle f, EntityHandle s, Range &seed_fe_range) const
base class for all interface classes